AUTOMATIC  TARGET  CLASSIFICATION  USING  HF  MULT I  FREQUENCY  1/2  ^ 
RADARS (U )  OHIO  STATE  UNIV  COLUMBUS  ELECTROSCIENCE  LAB 
J  CHEN  OCT  83  ESL-714190-3  N00014-82-K-0027 


UNCLASSIFIED 


F/G  17/9 


NL 


MICROCOPY  RESOLUTION  Ttbl  <  “ART 

NAHuNAl  "*  '-lANl.'AHl'-  > 


rogwrmn?rBTt?qT  n  w  v^»  \  :is  I  .'Vi.v  I  *»  t  -.-’.-■c’wi.ti  i~w  l~~»  u~v  U>'  ’jv  L"x  -"K  u>  •.  *  l-*  u 


The  Ohio  State  University 

AUTOMATIC  TARGET  CLASSIFICATION  USING 
HF  MULTI FREQUENCY  RADARS 

Jenshlun  Chen 


O) 

Tj- 

Mh 

CM 

CD 


< 

I 

D 

< 


er>- 


The  Ohio  State  University 

ElectroScience  Laboratory 

Department  of  Electrical  Engineering 
Columbus,  Ohio  43212 

Technical  Report  714190-3 
Contract  N00014-82-K-0037 
October  1983 


«.\0 


i 


Department  of  the  Navy 
Office  of  Naval  Research 
800  North  Quincy  Street 
Arlington,  Virginia  22217 


via  i  hi  nun  nr. 

—  -  '  ”  *  n  *  x,tvi£,|N  1  | 

Approved  tOJ  puhl.c  releaa^T  I 
-  Distribution  Unlilajte^  ^  J 


O 


NOTICES 


When  Government  drawings,  specifications,  or  other  data  are 
used  for  any  purpose  other  than  in  connection  with  a  definitely 
related  Government  procurement  operation,  the  United  States 
Government  thereby  incurs  no  responsibility  nor  any  obligation 
whatsoever,  and  the  fact  that  the  Government  may  have  formulated, 
furnished,  or  in  any  way  supplied  the  said  drawings,  specifications, 
or  other  data,  is  not  to  be  regarded  by  implication  or  otherwise  as 
in  any  manner  licensing  the  holder  or  any  other  person  or  corporation, 
or  conveying  any  rights  or  permission  to  manufacture,  use,  or  sell 
any  patented  invention  that  may  in  any  way  be  related  thereto. 


REPORT  DOCUMENTATION  »•  "EPORT  *• 

_ PAGE 

4.  TKI*  and  Subtitle  ~ 

AUTOMATIC  TARGET  CLASSIFICATION  USING  HF 
MULTIFREQUENCY  RADARS 


7.  Authors) 

Jenshiun  Chen 


9.  Performing  Organization  Nam#  and  Address 

The  Ohio  State  University  ElectroScience  Laboratory 
Department  of  Electrical  Engineering 
Columbus,  Ohio  43212 


11  s#trrpanf?fif?nHr,^*w^ntiaw" 

Office  of  Naval  Research 
800  North  Quincy  Street 
Arlington,  Virginia  22217 


is.  supptementfjn^io^rk  reported  in  this  report  was  also  used  as  a  di 
Ohio  State  University  Department  of  Electrical  Engineering 
for  the  degree  Doctor  of  Philosophy. 


IS.  Abatract  (Limit:  200  words) 


1.  Racifxent  t  Accttnon  No 


S.  Report  Data 

October  1983 


S.  Performing  Org.ru, »t. on  Rapt  No 

ESL  714190-3 


10.  Protect  /Task /Work  Unit  No 


11.  ContracffC)  or  Grant(G)  No 


N00014-82-K-0037 


13.  Type  of  Report  &  Period  Covered 

Technical  Report 


ssertation  submitted  to  the 
as  part  of  the  requirements 


v  ^~The  classification  of  radar  targets  such  as  aircraft  and  ships  using  lower 
resonance-region  radar  returns  has  been  of  significant  interest  in  recent  years.  The 
H.F.  band  is  in  the  resonance  region  of  such  targets.  The  probability  of 
mi sclassifi cation  depends  upon  the  post-processing  signal-to-noise  power  ratio. 

Current  techniques  for  measuring  and  processing  the  amplitude  and  phase  of  H.F. 
radar  returns  are  reviewed.  The  post-processing  SNR  is  determined  as  a  function  of 
coherent  observation  time.  Based  on  features  extracted  from  the  radar  returns,  a 
frequency-domain  Nearest-Neighbor  classification  algorithm  and  a  time-domain 
Cross-Correlation  classification  algorithm  are  designed.  A  set  of  radar  backscatter 
measurements  of  scale  model  ships  and  aircraft  Is  used  to  generate  a  data  base  for 
testing  the  classification  algorithms.  Statistical  techniques  are  applied  to  determine 
the  probability  of  misclassification  as  a  function  of  the  post-processing  SNR.  The 
probability  of  misclassification  for  a  reasonable  amount  of  coherent  observation  time  i 
determined.  .  ' 


17.  Document  An«l|r»li 


Descriptor* 


b.  Identifier* /Open  (ndeO  Tsrmt 


C.  COSAT  I  PI eld/Group 


IS.  Aveilebllity  Statement 


19.  Security  Cist*  (This  Report) 

_ Unclassified 

10.  Security  Clot.  (Thu  Rege) 

Unclassified 


21.  No  of 


(See  ANSl-X39.lt) 


OPTIONAL  PQRM  272  (A-77) 
(Formerly  NT  IS- 3  b) 
Department  of  Commerce 


Page 


C.  Classification  In  the  Frequency  Domain .  49 

D.  Classification  in  the  Time  Domain .  56 

VI.  GENERATION  OF  A  DATA  BASE .  61 

A.  Introduction .  61 

B.  Measurement  and  Calibration  .  62 

C.  Connection  of  Two  Frequency  Bands .  66 

D.  Interpretation  and  Smoothing .  72 

E.  Post-Processed  Radar  Returns .  81 

V.  TESTING  OF  CLASSIFICATION  ALGORITHMS .  83 

A.  Introduction .  83 

B.  Noise  Model  . .  85 

C.  Estimation  of  the  Probability  of 

Mi sclassifi cation  .  90 

D.  Classification  of  Ships .  95 

E.  Classification  of  Aircraft .  124 

VI.  CONCLUSIONS  AND  RECOMMENDATIONS  .  142 

REFERENCES .  148 

APPENDIX  A .  151 

APPENDIX  B .  179 


o-  ■  ’ p— r.  r? 


■-p  '  -• 


■3»  <■  ■  '■' 


t 


b 


LIST  OF  TABLES 


Table 


2.1.  Specifications  of  HF  Radars 


Page 

12 


2.2.  World  Meterological  Organization  Sea  State  [1] .  15 

2.3.  Additional  Transmission  Loss  Due  to  Sea  State 

(Barrick  [17]) .  17 

5.1.  A  List  of  Ships  Used  in  the  Classification 

Experiments . . .  97 

5.2.  A  List  of  Ships  Used  in  the  Classification 

Experiments .  101 

5.3.  A  Matrix  Showing  the  Distance  Between  any  Two 


Classes.  The  Entry  at  Row  1  and  Column  J  is  the 
N.N.  Distance  Between  Class  I  and  Class  J.  Only 
Amplitude  Data  are  Used  to  Calculate  the  Distance  ...  102 


5.4.  A  Matrix  Showing  the  Distance  Between  Any  Two  Classes. 

Only  the  Phase  Data  are  Used  to  Calculate  the  N.N. 

Distance . 103 

5.5.  A  Matrix  Showing  the  Distance  Between  Any  Two  Classes. 

Both  Amplitude  and  Phase  Data  are  Used . 104 

5.6.  A  Matrix  Showing  the  Cross-Correlation  Coefficient 

for  Any  Two  Classes . 106 


5.7.  Minimum  S/N  Required  (in  dB)  for  10  Percent  or  Less 
Classification  Errors.  For  Each  Algorithm  (e.g., 

A&W)  There  are  Four  S/N  Values  (e.g.,  19,  17,  13, 

15)  Corresponding  to  Four  Different  Conditions 

Regarding  Aspect  Angles  and  the  Number  of  Frequencies 

Used  N . 119 


v 


v^.  V  ^  \  -V 

.  %  i W WLAT -*■'  -  a  -V  — »  -  »  -*■> 


»A,k  * < 


m*'M 


LIST  OF  FIGURES 


1.1.  Block  diagram  of  a  typical  target  classification 

system.  The  catalog  contains  radar  returns  of 
some  preselected  targets  (A,  6).  The  output  of  the 
signal  processor  is  a  set  of  amplitude  and  phase 
returns  of  the  unknown  target . 

2.1.  Family  of  transmission  curves  parametric  in 
frequency  for  a  fixed  distance  of  2,000  km 
superimposed  on  a  vertical  ionogram  (Davis  [5]). 

The  abscissa  is  the  vertical  transmission  frequency, 
and  the  parameter  is  the  oblique  transmission 
frequency . 

2.2.  Ray  paths  corresponding  to  intersections  of 
transmission  curves  and  the  ionogram  (Davis  [5])  .  . 

2.3.  Additional  transmission  loss  due  to  sea  state  at 

10  MHz  (Barrick  [17]) . 

2.4.  A  diffraction  grating  model  for  sea  surface 

waves . 

2.5.  A  Doppler  spectrum  is  shown  in  which  a  ship  target 

return  is  separated  from  the  clutter.  The 
reference  here  is  the  larger  first  order  sea 
scatter.  The  Doppler  frequency  here  is 
normalized  to  the  1st  order  sea  scatter  Doppler 
frequency . 

2.6.  The  area  illuminated  by  a  radar  is  shown  as  the 
shaded  area.  Note  here  that 

(f>  =  grazing  ange 
t  =  pulse  width 
C  =  speed  of  light 
R  =  wavepath  length 

3n  =  antenna  azimuth  beam  width . 


2.7.  Block  diagram  of  a  coherent  radar  system . 

3.1.  Basic  process  of  a  pattern  classification 

system . 

4.1.  Diagram  of  the  radar  measurement  facility  . 

4.2.  Simulation  of  the  reflection  from  the  sea  surface. 

A  replica  is  attached  to  each  ship  model . 

4.3.  Connection  of  two  frequency  bands  . 

4.4.  Determination  of  the  relative  position  offset  in 

the  time  domain  . 

4.5.  Interpolation  and  smoothing  of  a  noisy  waveform 

using  a  Hamming  Window . 

4.6.  Interpretation  of  the  interpolation  in  the  time 

domain.  The  Hemming  window  becomes  a  low-pass 
filter  in  the  time  domain  . 

4.7.  An  example  of  the  interpolation.  There  are  6  data 
points  located  at  f  =  1,  3,  3.5,  4,  4.5,  and  5. 

The  underlying  function  is  c(f)  =  f.  It  is  desired 
to  estimate  the  data  value  at  f  =  3  . 

4.8.  Comparison  between  measured  and  theoretical  returns 

of  a  3.5  inch  sphere  at  f  *  2-4GHZ . 

5.1.  A  flowchart  showing  the  experimental  process  for 

testing  the  classification  algorithms  . 

5.2.  The  distribution  of  the  noise  on  an  I-Q  plane  .  .  . 

5.3.  Change  of  wavepath  length  Ap  for  a  change  in 
aspect  angles  by  15°.  Note  Ap  »  2a«.  =  2. 

(D/2  -  D/2  cos as)  . 

5.4.  Plots  of  error  probabilities  with  the  number  of 

frequencies  N  as  a  parameters.  Only  amplitude  is 
used.  Aspect  angles  are  assumed  known.  A  90% 
confidence  interval  for  a  classification  error  of 
30%  is  approximately  ±5%  . 


i>  >  .•  yftj^ 


Figure 


Page 


^  * 


ft 


0 


5.5.  Plots  of  error  probabilities  with  the  number  of 
frequencies  N  as  a  parameters.  Only  phase  is  used. 

Aspect  angles  are  assumed  known.  A  90%  confidence 
interval  for  a  classification  error  of  30%  is 
approximately  ±5% .  110 

5.6.  Plots  of  error  probabilities  with  the  number  of 
frequencies  N  as  a  parameters.  Both  amplitude  and 
phase  are  used.  Aspect  angles  are  known  a  priori. 

A  90%  confidence  interval  for  a  classification 

error  of  30%  is  approximately  t  5% .  Ill 

5.7.  Error  probabilities  for  N=2,4,8  frequencies.  Time- 
domain  algorithm  is  used,  and  the  aspect  angles  are 
assumed  a  priori.  A  90%  confidence  interval  for  a 


classification  error  of  30%  is  approximately  ±  5%.  .  .  112 

5.8.  Error  probabilities  for  four  algorithms.  Number  of 
frequency  =  8.  Aspect  angles  are  known  a  priori. 

A  90%  confidence  interval  for  a  classification  error 

of  30%  is  approximately  t  5% .  113 

5.9.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  4.  Known  aspect  angles.  A  90% 
confidence  interval  for  a  classification  error  of 

30%  is  approximately  t  5% .  114 

5.10.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  2.  Aspect  angles  are  known  a  priori. 

A  90%  confidence  interval  for  a  classification  error 

of  30%  is  approximately  t  5% .  115 

5.11.  Error  probabilities  for  relative  amplitudes.  Numbers 
of  frequencies  =  8.  Aspect  angles  are  known.  A 

90%  confidence  interval  for  a  classification  error 

of  30%  is  approximately  t  5% .  117 


5.12.  Comparison  of  error  probability  between  the  cases 
where  aspect  angles  are  either  known  or  unknown. 

Only  amplitude  is  used.  N  =  8  frequencies.  A 

90%  confidence  interval  for  a  classification  error 

of  30%  is  approximately  t  5% .  120 


i  x 


§ 


ii. 


Figure 


Page 


5.13.  Comparison  of  error  probalities  between  the  cases 
for  known  aspect  angles  and  unknown  aspect  angles. 

Only  phase  is  used.  N  =  8  frequencies.  A  90% 
confidence  interval  for  a  classification  error  of 

30%  is  approximately  ±5% .  121 

5.14.  Error  probabilities  for  "known  aspect  angles"  and 
"unknown  aspect  angles."  Both  amplitude  and  phase 
are  used.  Number  of  frequencies  =8.  A  90% 
confidence  interval  for  a  classification  error  of 


30%  is  approximately  t  5% .  122 

5.15.  Error  probabilities  for  the  cases  where  the  aspect 
angles  are  either  known  or  unknown.  Time-domain 
algorithm  is  used.  Number  of  frequencies  =  8. 

A  90%  confidence  interval  for  a  classification 

error  of  30%  is  approximately  t  5% .  123 

5.16.  Error  probabilities  for  N  =  2,4,8  freququencies.  Only 

amplitude  data  are  used.  The  aspect  angles  are  assumed 
known.  A  90%  confidence  interval  for  30%  classification 
error  is  approximately  t  5.6% .  127 

5.17.  Error  probabilities  for  N  *  2,4,8  frequencies. 

Only  phase  data  are  used.  The  aspect  angles  are 
assumed  known.  A  90%  confidence  interval  for  30% 
classification  error  is  approximately  t  5.6%  .  128 

5.18.  Error  probabilities  for  N  =  2,4,8  frequencies.  Both 
amplitude  and  phase  data  are  used.  The  aspect  angles 
are  assumed  known.  A  90%  confidence  interval  for 

30%  classification  error  is  approximately  t  5.6%  .  .  .  129 


5.19.  Error  probabilities  for  N  =  2,4,8  frequencies.  The 

time-domain  algorithm  is  used.  The  aspect  angles  are 
assumed  known.  A  90%  confidence  interval  for  30% 
classification  error  is  approximately  -  5.6%  .  130 

5.20.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  8,  and  the  aspect  angles  are  known. 

A  90%  confidence  interval  for  30%  classification 

error  is  approximately  t  5.6% .  134 


x 


Fi gure 


Page 


i 


m 

t  ■ 

r, 

*l 


r * 


* 


r'- 


* 


5.21.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  8,  and  the  aspect  angles  are  known. 

A  90%  confidence  interval  for  30%  classification 

error  is  approximately  ±5.6% .  135 

5.22.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  2,  and  the  aspect  angles  are  known. 

A  90%  confidence  interval  for  30%  classification 

error  is  approximately  ±  5.6% .  136 

5.23.  Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  8,  and  the  aspect  angles  are  known. 

A  90%  confidence  interval  for  30%  classification  error 
is  approximately  t  5.6% .  137 

5.24.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  Only 
amplitude  data  are  used,  and  number  of  frequencies  =  8. 

A  90%  confidence  interval  for  30%  classification  error 

is  approximately  ±5.6%  .  138 


5.25.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  Only  phase 
data  are  used,  and  number  of  frequencies  =  8. 

A  90%  confidence  interval  for  30%  classification 

error  is  approximately  ±  5.6% .  139 

5.26.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  Both 
amplitude  and  phase  data  are  used,  and  number  of 
frequencies  =8.  A  90%  confidence  interval  for 

30%  classification  error  is  approximately  ±  5.6%  .  .  .  140 

5.27.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  The  time- 
domain  algorithm  are  used,  and  number  of  frequencies  =  8. 

A  90%  confidence  interval  for  30%  classification  error 

is  approximately  ±  5 . 6% .  141 


A.  1 .  Ship  A  at  0  degrees .  152 

A. 2.  Ship  A  at  90  degrees .  153 

A. 3.  Ship  A  at  180  degrees .  154 

A. 4.  Ship  B  at  0  degrees . 155 


xi 


I 


Figure  Page 

A. 5.  Ship  B  at  90  degrees .  156 

A. 6.  Ship  B  at  180  degrees .  157 

A. 7.  Ship  C  at  0  degrees .  158 

A. 8.  Ship  C  at  90  degrees .  159 

A. 9.  Ship  C  at  180  degrees .  160 

A. 10.  Ship  D  at  0  degrees .  161 

A.  11.  Ship  0  at  90  degrees .  162 

A.  12.  Ship  0  at  180  degrees .  163 

A.  13.  Ship  E  at  0  degrees . 164 

A.  14.  Ship  E  at  90  degrees .  165 

A.  15.  Ship  E  at  180  degrees .  166 

A.  16.  Aircraft  F  at  0  degrees .  167 

A.  17.  Aircraft  F  at  90  degrees .  168 

A.  18.  Aircraft  F  at  180  degrees .  169 

A.  19.  Aircraft  G  at  0  degrees .  170 

A. 20.  Aircraft  G  at  90  degrees .  171 

A. 21.  Aircraft  G  at  180  degrees .  172 

A. 22.  Aircraft  H  at  0  degrees .  173 

A. 23.  Aircraft  H  at  90  degrees .  174 

A. 24.  Aircraft  H  at  180  degrees .  175 

A. 25.  Aircraft  I  at  0  degrees .  176 

A. 26.  Aircraft  I  at  90  degrees .  177 

A. 27.  Aircraft  I  at  180  degrees .  178 


*  * 

;a 


.-•i 


xii 


CHAPTER  I 


INTRODUCTION 

A.  GENERAL  NATURE  OF  THE  PROBLEM 

Radars  operating  at  HF  (3  to  30  MHz)  have  a  unique  capability  of 
allowing  propagation  over  long  distance.  Skolnik  and  Headrick  [1,  2] 
indicated  that  over-the-horizon  (OTH)  HF  radars  using  sky-wave 
propagation  via  refraction  by  the  ionoshere  are  capable  of  detecting 
such  targets  as  aircraft, ships,  and  missiles  at  distances  beyond  4,000  km. 
The  wavelengths  of  the  HF  band  range  from  10  to  100  m  and  are 
comparable  to  the  dimensions  of  such  targets.  This  means  that  the  HF 
band  is  in  the  resonance  region  of  such  targets.  Kennaugh  and  Moffatt 
[3,  4]  found  that  low-frequency  radar  returns  in  the  resonance  region 
carried  essential  information  regarding  the  overall  dimensions  and 
approximate  shapes  of  the  targets.  Based  on  their  work  much  research 
effort  has  been  directed  toward  the  identification  of  targets 
using  resonance-region  returns  [4  -  13].  Ksienski  and  Lin  [5,  8,  10] 
showed  that  reliable  classification  of  aircraft  can  be  achieved  by 
processing  resonance-region  mul ti-frequency  radar  returns,  and 
Walton  and  Young  [13]  demonstrated  similar  results  for  ships. 
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Fig.  1.1  shows  a  diagram  of  an  HF  target  classification  system. 
Possible  targets  include  surface  ships,  aircraft,  and  missiles.  The 
radar  is  capable  of  operating  at  several  frequencies  in  the  HF  band. 

Two  mechanisms  of  wave  propagation  are  possible:  sky  wave  and  surface 
wave.  The  range  of  detection  may  vary  depending  on  the  propagation 
mechanism,  frequency,  weather,  and  other  factors.  The 
multi -frequency  radar  returns  detected  by  the  receiver  are  processed, 
and  the  target  classifier  then  attempts  to  classify  the  detected  target 
or  targets  based  on  the  processed  radar  returns.  The  input  to  the 
target  classifier  is  a  set  of  amplitude  and  phase  measurements  of  the 
target  returns  at  selected  frequencies.  Conceptually  the  target 
classifier  first  extracts  features  from  the  radar  returns  then 
'compares'  the  features  with  those  in  the  catalog  or  data  base,  which  is 
a  set  of  ideal  radar  returns  of  a  number  of  pre-selected  targets.  The 
target  classifier  can  either  identify  an  unknown  target  as  one  of  the 
cataloged  targets  or  claim  that  it  is  not  listed  in  the  catalog.  The 
probability  of  mi sclassifi cation  depends  on  the  post-processing 
signal-to-noise  ratio  of  the  radar  returns,  and  the  classification 
schemes  used. 

B.  RESEARCH  OBJECTIVES 

Previous  research  work  (e.g.  ref.  4  -  13)  usually  assumed  that 
amplitude  and/or  intrinsic  phase  were  available.  For  HF  radars  the 
amplitude  returns  may  vary  widely  as  a  function  of  the  propagation  lass 


Figure  1.1.  Block  diagram  of  a  typical  target  classification 
system.  The  catalog  contains  radar  returns  of 

A  A 

some  preselected  targets  (A,  0).  The  output  of  the 
signal  processor  is  a  set  of  amplitude  and  phase 
returns  of  the  unknown  target. 
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antenna  gains,  and  other  factors.  The  phase  returns  are  a 
function  of  wave  path  length  and  frequency,  and 

can  vary  over  the  entire  branch  (e.g  0°  to  360°).  Therefore  such  radar 
returns  must  be  processed  extensively  before  reliable  features  can  be 
extracted  to  classify  the  unknown  targets.  The  probability  of 
correctly  classifying  the  unknown  targets  depends  on  the 
accuracy  with  which  the  radar  returns  are  measured  and  the 
classification  algorithms  employed.  The  major  purpose  of  this  research 
is  to 

1.  Study  current  techniques  for  measuring  and  processing 
the  target  returns. 

2.  Determine  the  approximate  values  of  the  post-processing 
signal-to-noise  power  ratio  as  a  function  of  coherent 
observation  time. 

3.  Develop  target  classification  algorithms  which  utilize 
the  processed  amplitude  and  relative  phase  returns. 

4.  Determine  the  approximate  values  of  the  probability  of 
mi sclassifi cation  as  a  function  of  post-processing 
signal-to-noise  power  ratio. 

5.  Determine  the  probability  of  misclassification 
achievable  for  a  reasonable  amount  of  coherent 


observation  time. 


C.  ORGANIZATION  OF  THIS  DISSERTATION 


Chapter  II  will  discuss  two  types  of  HF  radars,  namely,  sky-wave 
and  surface-wave  radars  [1,  2,  14].  The  ionospheric  and  surface 
propagation  effects  will  be  briefly  discussed  [15,  16,  17].  Current 
techniques  for  measuring  the  amplitude  returns  [18]  and  relative  phase 
returns  [19]  will  also  be  reviewed.  An  estimation  of  the 
post-proscessing  signal-to-noise  power  ratio  as  a  function  of  total 
coherent  observation  time  will  be  derived  [20]. 

Chapter  III  will  begin  with  a  discussion  of  the  basic 
classification  process,  and  a  brief  discussion  of  Bayes  decision  theory 
[21,  22].  Based  on  the  results  of  Chapter  II,  two  non-parametric 
classification  algorithms  will  be  presented.  The  first  algorithm  will 
make  use  of  the  Nearest-Neighbor  decision  rule  to  classify  targets  in 
the  frequency  domain  [5,  8,  10,  13].  The  Nearest-Neighbor  rule  had  been 
studied  extensively  in  the  literature  and  was  proved  to  be  very 
effective  in  identifying  targets  [5,  8,  10,  13].  The  second  will  use 
cross-correlation  processing  to  classify  targets  in  the  time  domain  [3, 
9].  The  correlation-processing  is  basically  Matched-filtering,  and  has 
been  used  extensively  in  signal  detection  but  has  seldom  been  applied  to 
target  identification.  Both  algorithms  require  a  catalog  as  mentioned 
above. 

Chapter  IV  will  discuss  techniques  for  generating  a  data  base  [13]. 
Techniques  for  processing  measured  radar  returns  will  be  presented. 

Radar  measurements  of  scale  model  targets,  which  were  conducted  in  an 


anechoic  chamber,  will  be  discussed.  The  processing  of  the  measured 
data  and  the  generation  of  an  experimental  database  using  the  processed 
data  will  be  discussed. 

Chapter  V  will  discuss  testing  of  the  classification  algorithms 
using  the  experimental  database.  The  targets  selected  include  ships  and 
aircraft.  This  chapter  will  estimate  the  probability  of 
mi sclassifi cation  as  a  function  of  signal-to-noise  ratio.  All  targets 
will  be  assumed  to  be  listed  in  the  database  [23].  The  design  of  a 
noise  model  and  the  scheme  for  estimating  the  classification  probability 
using  a  finite  number  of  samples  will  also  be  discussed  [24].  The 
probability  of  mi sclassifi cation  as  a  function  of  coherent  observation 
time  will  be  determined. 

Chapter  VI  will  present  a  summary  of  the  significance  of  the  work, 
conclusions,  and  recommendations  for  future  work. 


;v> 


CHAPTER  II 


MEASUREMENT  OF  H.F.  RADAR  RETURNS 


A.  INTRODUCTION 


As  discussed  in  Chapter  I,  our  approach  to  target  classification  is 
to  compare  features  extracted  from  radar  returns  of  a  target  with  those 
in  the  data  base.  It  is  evident  that  the  probability  of  correctly 
classifying  the  targets  is  closely  related  to  the  accuracy  with  which  we 
can  measure  the  radar  returns.  Intuitively  we  know  that  the  higher  the 
signal -to-noise  ratio,  the  lower  the  error.  The  use  of  HF  frequencies 
for  long  distance  detection  further  complicates  the  problem  because  of 
the  large  amount  of  noise  that  will  accompany  the  radar  returns. 
Therefore  the  design  of  identification  algorithms  must  be  facilitated  by 
a  knowledge  of  the  characteristics  of  the  post-processed  H.F.  radar 
returns.  The  purpose  of  this  chapter  is  to  examine,  in  general,  how 
well  we  can  estimate  the  H.F.  radar  returns.  The  post-processing 
signal -to-noise  ratio  will  be  estimated  as  a  function  of  coherent 


observation  time. 


This  will  enable  us  to  acquire  an  estimation  of  the  probability 
of  misclassification  as  a  function  of  coherent  observation  time,  and 
will  give  us  an  idea  about  the  error  probability  we  will  likely  have  in 
practice. 

Section  B  will  discuss  the  characteristics  of  two  types  of  HF 
radars,  i.e.,  sky-wave  and  surface-wave  HF  radars,  and  a  brief 
discussion  on  the  propagation  effects. 

Section  C  will  address  the  problem  of  measuring  the  radar  cross 
section  (RCS)  of  a  target.  It  will  examine  and  discuss  the  current 
techniques  for  calibrating  the  target  cross  section  using  the  sea 
scatter  as  a  reference.  The  post-processing  signal-to-noise  (SNR)  power 
ratio  will  be  evaluated  statistically. 

Section  D  will  discuss  measurement  of  the  phase  of  radar  returns 
using  coherent  radars,  and  determine  the  constituents  of  the  coherent 
phase  returns  at  the  output  of  the  receiver. 

B.  SOME  CHARACTERISTICS  OF  HF  RADARS 

This  section  will  first  discuss  sky-wave  OTH  H.F.  radars  [1,  2] 
and,  then,  discuss  surface-wave  OTH  H.F.  radars  [14]  and  a  comparison 
between  the  two  types  of  radars. 

1 .  H.F.  Sky-Wave  Radars 

H.F.  sky-wave  radars  make  use  of  ionospheric  refraction  for  the 
detection  of  targets  at  over-the-horizon  distances.  Ionospheric 
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conditions  may  vary  with  weather,  season,  time  of  day,  solar  activity, 
etc.,  and  the  optimum  operating  frequencies  may  change  from  time  to  time 
over  a  wide  portion  of  the  H.F.  band. 

There  may  be  more  than  one  refracting  region  in  the  ionosphere. 

The  highest  region,  i.e.,  the  F2  layer,  is  at  altitudes  of  from  230  to 
400  km.  A  single  refraction  from  such  region  allows  radar  ranges  to  be 
extended  to  almost  4000  km.  Consider  a  figure  from  [15],  as  shown  in 
Fig.  2.1.  The  figure  is  a  diagram  of  transmission  curves  with  the 
oblique  transmission  frequency  as  a  parameter  for  a  fixed  distance  of 
of  2,000  km  superimposed  on  a  vertical  ionogram.  Here  the  transmission 
curves  describe  the  relationship  between  the  oblique  transmission  and 
equivalent  vertical  transmission  frequencies.  Note  that  the  equivalent 
vertical  transmission  frequency  is  the  frequency  of  the 
vertically-transmitted  waves  which  reach  the  same  virtual  height  as  the 
obliquely-transmitted  waves  do.  In  this  case  a  plane  reflector  is 
used  to  model  the  ionosphere,  and  the  transmission  curves  are  given  by 
the  'Secant  Law'  [15],  i.e., 

^ob  =  fv  sec  <t.0 


where  h'  =  equivalent  virtual  height 
D  =  distance 

f0b  =  oblique  transmission  frequency 


VIRTUAL  HEIGHT  ,  Kr 


.  Family  of  transmission  curves  parametric  in  frequency  for 
fixed  distance  of  2,000  km  superimposed  on  a  vertical 
ionogram  {Davis  [15]).  The  abscissa  is  the  vertical 
transmission  frequency,  and  the  parameter  is  the  oblique 
transmission  frequency. 


fv  =  vertical  transmission  frequency 
=  elevation  angle. 

The  intersection  of  a  transmission  curve,  for  a  qiven  oblique 
transmission  frequency  f0b,  with  the  vertical  ionogram  gives  the  virtual 
heights  of  reflection  of  the  transmitted  curves.  In  Fig.  2.1  the  three 
sections  of  the  vertical  ionogram  would  correspond  to  reflection  from 
the  E,  F i  and  F2  layers  of  the  ionosphere.  A  frequency  of  14  MHz, 
for  example,  for  propagation  over  a  distance  of  2,000  km,  could  travel 
by  several  different  paths,  two  of  which  would  be  reflected  from  the 
F2  layer  at  heights  corresponding  to  the  intersections  marked  a  and  a' 
on  the  vertical  ionogram.  The  corresponding  ray  paths  are  shown  in  Fig. 
2.2  [15].  For  values  of  fv  above  20  MHz  there  is  no  intersection, 
and,  hence,  no  level  of  reflection.  The  distance  for  which  a  given 
frequency  is  the  maximum  frequency  is  called  the  skip  distance. 

Within  this  distance  no  signals  are  received  from  the  ionosphere. 

Waves  with  different  frequencies  may  have  different  phase  delays 
because  they  may  propagate  along  different  paths.  The  extraction  of 
features  from  mul ti -frequency  phase  returns  must  therefore  contend 
with  this  frequency-dispersion  phenomenon. 

In  general ,  it  is  difficult  to  predict  the  received  signal  level 
because  it  may  require  a  detailed  knowledge  of  the  ionosphere  along  the 
ray  path  together  with  elaborate  ray-tracing  [15].  Navon  and  Cory  [16] 
showed  that  the  received  signal  median  for  a  one-way  communication  link 
of  about  2000  km  can  be  predicted  with  an  accuracy  of  about  5  dB  which 
is  unfortunately ,  too  high  to  be  useful  for  target  identification. 
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TABLE  2.1 


SPECIFICATIONS  OF  HF  RADARS 


SKY-WAVE  RADAR 

SURFACE -WAVE  RADAR 

TRANSMITTED  POWER 

SEVERAL  HUNDRED  kW 

HIGHER 

RANGE 

1000  to  4000  km 

about  300  km 

Antenna  gain 

20  to  30  dB 

same 

Antenna  horizontal 
length 

about  1  km 

about  1  km 

Waveform 

can  be  pulsed 
sinusoidal  waves 

can  be  pulsed 
sinusoidal  waves 

Bandwidth 

5  to  100  kHz 

can  be  wider 

Pulse  width 

200  us  to  10  us 

can  be  narrower 

Range  resolution 

2  to  40  km 

can  be  smaller 

Doppler  resolution 

about  0.1  Hz 

can  be  less 
than  0.1  Hz 

PRF 

about  30  Hz 

about  300  Hz 

Targets 

aircraft,  ships, 
mi ssiles 

same 

Antenna  Azimuth 
beamwi dth 

about  1° 

about  1° 

Fig.  2. 


l*I.AfU:  lONOSHIII-'.KF. 


2.  Ray  paths  corresponding  to  intersections  of  transmission 


1 


5 


curves  and  the  ionogram  (Davis  [15]). 


The  range  of  operation  for  sky-wave  radar  can  be  extended  up  to 
4,000  km  and  beyond.  However  the  ranqe  of  operation  is  generally 
limited  to  a  minimum  of  about  1,000  km  due  to  the  ionospheric 
skip-distance  effect  [15].  A  specification  of  a  typical  HF  sky-wave 
radars  was  qiven  by  Skolnik  and  Headrick  [1,  2],  and  is  shown  in  Table 
2.1. 

2 .  H.F.  Surface-Wave  Radars 

Surface-wave  radars  rely  on  surface  or  ground  wave  propagation 
around  the  curvature  of  the  earth.  Because  of  high  losses,  propagation 
over  long  distances  can  best  be  achieved  over  the  ocean  surface  [141. 

The  range  of  operation  can  vary  from  tens  of  kilometers  up  to  about  500 
kilometers.  It  covers  the  range-gap  between  line-of-siqht  radars  and 
sky-wave  radars.  The  specification  for  a  typical  surface  wave  radar  is 
quite  different  from  that  of  a  sky-wave  radar,  as  can  be  seen  from  Table 
2.1.  This  will  be  briefly  explained  below. 

Barrick  [17]  indicated  that  the  sea  state  can  significantly  affect 
the  ground-wave  loss  for  propagation  across  the  ocean  at  HF.  The  sea 
state  is  a  measure  of  the  wave  height,  as  shown  in  Table  2.2  [11. 

Figure  2.3  is  a  plot  of  the  additional  transmission  loss  due  to  sea 
state  at  10  MHz  [17].  It  can  be  seen  that  the  transmission  loss 
increases  rapidly  as  the  ground  distance  increases.  Table  2.3  is  a  list 
of  additional  transmission  loss  due  to  sea  state  for  a  fixed  ground 
distance  of  463  km  (250  nmi).  It  can  be  seen  that,  for  sea  states  of  4 
and  5,  the  two-way  losses  at  a  frequency  of  10  MHz  can  be  as  high  as  23 


TABLE  2.2 


World  Meteorological  Organization  Sea  State  [1]. 


Wave  height 


Sea  state 

Feel 

Meters 

Descriptive  term 

0 

0 

0 

Calm,  glassy 

1 

o-i 

0-0  1 

Calm,  rippled 

2 

i-H 

01-0  5 

Smooth,  wavelets 

3 

2-4 

0  6-1.2 

Slight 

4 

4-8 

12-2  4 

Moderate 

5 

8-13 

2  4-4  0 

Rough 

6 

13-20 

4.0-60 

Very  rough 

7 

20-30 

6.0-90 

High 

8 

30-45 

90-14 

Very  high 

9 

over  45 

over  14 

Phenomenal 

TABLE  2.3 


FREQUENCY 

(MHz) 


SEA  STATE  =  4 


SEA  STATE 


dB  and  31  dB,  respectively,  above  the  free-space  loss.  It  can  also  be 
seen  that  the  loss  is  less  for  low  frequencies.  Because  of  high 
propagation  loss,  surface-wave  radars  have  a  much  shorter  range  of 
detection  than  sky-wave  radars. 

In  most  cases  the  incident  elevation  angle  of  surface  waves  on 
targets  is  equal  to  or  slightly  greater  than  zero,  whereas  the  incident 
elevation  angles  of  sky  waves  may  vary  widely,  depending  on  target 
range.  The  latter  may  require  a  much  larger  data  base  containing 
backscatter  data  for  a  set  of  incident  elevation  angles,  resulting  in 
tremendous  increase  in  complexity,  ambiguity  in  look  angles,  and  costs. 

Finally  in  this  section  we  would  like  to  note  that  H.F.  radars  must 
employ  doppler  processing  to  seperate  desired  moving  targets  from 
surface  clutter,  because  the  clutter  returns  are  usually  much  larger 
than  the  target  returns  [1],  The  detection  of  low-speed  targets  such  as 
ships  requires  more  sophisticated  processing  and  longer  processing  time 
to  provide  sufficient  resolution  in  doppler  frequency.  Since  a 
processing  time  of  T  seconds  corresponds  to  a  doppler  frequency 
resolution  of  1/T  HZ,  we  can,  ideally,  obtain  a  very  high  doppler 
resolution  if  we  increase  the  processing  time.  Trizna  [5],  however, 
indicated  that  there  is  a  practical  limit  on  the  minimum  spectral  width 
which  can  be  achieved. 

He  indicated  that  the  most  severe  resolution  limitation  is  provided 
by  the  ionosphere.  Propagation  in  the  ionosphere  results  in  phase 
changes  which  may  not  be  uniform  across  the  phase  fronts  of  the  radar 
waves.  The  ionosphere  condition  may  change  with  time,  resulting  in 


Eq.  (2.6)  shows  that  the  amplitude  return  is  a  function  of  many 
parameters,  some  of  which  may  be  difficult  to  estimate.  For  example  the 
propagation  loss  Lp  for  either  sky-wave  or  surface  wave  radars  is 
difficult  to  estimate  since,  as  discussed  in  section  B,  it  may  vary  as 
the  ionospheric  condition  or  sea  state  changes.  The  amplitude  as 
calculated  by  Eq.  (2.6)  cannot  be  used  directly  for  target 
classification.  One  way  to  reduce  the  number  of  parameters  is  to  use  as 
a  reference  the  signal  return  from  another  object  which  is  close  to  the 
target,  and  whose  radar  cross  section  can  be  estimated  independently. 
Suppose  the  reference  signal  power  at  the  receiver  output  is  given  as 
c  ref,.  \3n4.  2, 

aref  =  p  ( 4tt )  R  Lp  Ls  (2.7) 

PTGTGRX  Ga 

where  oref  is  the  radar  cross  section  of  the  reference  object.  If  the 
reference  objects  are  close  enough  to  the  target,  the  propagation  loss 
Lp  and  the  wave  path  length  R  in  Eq.  (2.7)  should  be  close  to  those  in 
Eq.  (2.5).  Other  parameters  Py,  Gy,  Gr,  Ga,  and  Ls  in  both  equations 
should  be  similar  because  the  target  and  the  reference  are  illuminated 
by  the  same  radar  at  the  same  time.  Taking  the  ratio  of  Eq.  (2.5)  and 
Eq.  (2.7),  we  have 

pr  =,  o 
preT  ^ref 

or 


a  13 


(  Pr  )  o^e f 


(2.8) 


Apparently  Eq.  (2.8)  above  is  much  simpler  and  more  reliable  than  Eq. 
(2.5).  Note  that  Eq.  (2.8)  is  a  calibration  equation  for  a  particular 
frequency.  For  multi -frequency  returns  the  returns  of  each  frequency 
must  be  calibrated  individually.  The  parameters  pr  and  pre^  can  be 
measured  using  Doppler  processing,  and  are^  can  be  predicted,  as 
discussed  below. 

Trizna  [18]  suggested  that  the  sea  scatter  can  be  used  as  a 
reference  for  calibrating  the  target  cross  section.  The  sea  may  be 
regarded  as  composed  of  a  large  number  of  individual  wave  trains,  each 
with  a  different  wavelength  and  amplitude  and  travelling  in  different 
directions  [1].  Assuming  that  the  sea  acts  as  a  simple  diffraction 
grating,  Crombie  [25]  showed  that,  for  radar  waves  incident  on  the  sea 
surface  with  an  incident  angle  of  <j>,  the  radar  responds  mainly  to  two 
resonant  sea  wave  trains  which  are  travelling  radially  toward  and  away 
from  the  radar,  as  shown  in  Fig.  2.4.  From  the  figure,  it  is  evident 
that  the  scattered  signals  are  in-phase  when 

2X<d  cos<|>  =  nx  ;  n  =  1 ,  2,  3  .  .  . 


where 


Xu  =  water  wave  length 
X  =  radar  wave  length 
4>  =  wave  incident  angle 


The  sea  scatter  is  called  first-order  for  n  =  1,  second-order  for  n  = 


Fig.  2.4.  A  diffraction  grating  model  for  sea  surface  waves 


and  so  forth.  The  velocity  of  a  water  wave  which  is  driven  by  gravity 
[26]  is  given  by 


where  g  =  acceleration  due  to  gravity.  Therefore  the  Doppler 
frequencies  of  the  water  waves  are  given  by 

fd  =  2v  =  ±  /  ng  (2.9) 

X  \J  ttXC0S4> 

After  Crombie  identified  the  mechanism  of  sea  scatter,  considerable 
research  has  been  conducted  on  the  analysis  of  sea  scatter  [18,  20, 

26,  27,  28,  29,  30].  Now  let's  discuss  the  calibration  of  target  cross 
section  using  sea  scatter,  as  suggested  by  Trizna  [18]. 

For  convenience  a  figure  from  [18]  is  shown  in  Fig.  2.5,  which  is  a 
radar  spectrum  of  sea  echo  collected  with  »  17-second  integration  time 
and  a  21.8-MHZ  radar  frequency,  at  a  range  of  about  1600  KM.  Note  that 
the  frequency  in  the  figure  is  normalized  to  the  first-order  sea 
scatter.  The  area  illuminated  by  the  radar  for  a  pulse  length  of  x  and 
azimuth  beam  width  63  is  given  by  [1]  {see  Fig.  2.6) 

Ac  =  R0B  (£i)  sec$ 

2 

where 

Ac  =  illuminated  area 
R  =  wave  path  length 
=  horizontal  antenna  beam  width 


DOPPLER  FREQUENCY  IV 


Fig.  2.5.  A  Doppler  spectrum  is  shown  in  which  a  ship  target 

return  is  separated  from  the  clutter.  The  reference  here  is 
the  larger  first  order  sea  scatter.  The  Doppler  frequency 
here  is  normalized  to  the  1st  order  sea  scatter  Doppler 
frequency. 


4>  =  elevation  angle 
x  =  pulse  width 

The  illuminated  area  (or  scattering  cell)  in  this  case  was  of  the  order 
of  90  dB-square  meters.  As  shown  in  the  Doppler  spectrum,  the  ship  return 
is  separated  from  the  sea  scatter,  and  the  target  return  and  the  sea 
scatter  can  be  readily  obtained.  It  can  be  seen  that  the  value  of  the 
cross  section  of  the  target  is  22  dB  below  the  greater  first-order  sea 
scatter  located  at  f  =  -1.  Based  on  a  second-order  model,  Trizna  [18] 
theoretically  predicted  that  the  magnitude  of  the  sea  scatter  in  this 
case  is  51  dB-square-meters.  Therefore  he  found  that  the  target  cross 
section  is  29.2  dB-square-meters  or  832  square  meters. 

Trizna  indicated  that  the  calibrated  target  cross  section  was 
sufficiently  accurate  to  be  applied  to  most  conditions  appropriate  to  HF 
propagation,  in  particular  the  high  frequency  portion  of  HF  band.  As 
the  techniques  are  being  updated  and  refined,  the  accuracy  of  measuring 
the  target  RCS  will  gradually  improve. 

The  radar  cross  section  of  the  sea  is  a  function  of  wave  height, 
current  direction,  wind  speed,  etc.,  and  will  vary  as  weather  and  sea 
state  change  [17,  18].  Barrick  [20]  experimentally  found  that  the 
spectral  component  of  the  first-order  sea  scatter  can  be  approximated  by 
a  Gaussian  distribution.  This  enables  us  to  analyze  the  error  embedded 
in  the  calibrated  target  cross  section  a  in  Eq.  (2.8),  as  discussed 


In  order  to  increase  the  accuracy  of  the  estimated  target  cross 
section  we  can  increase  the  number  of  observations  and  then  calculate 


the  average  of  the  received  radar  returns.  Let  p  be  the  sample  average 
defined  as 


M 

p  =  1  I  Pi  (2.10) 

R  i=l 

where  P-j  is  the  received  power  for  the  ith  observation.  From  Eq.  (2.8) 
and  Eq.  (2.10)  we  can  estimate  the  target  cross  section  a  as 


where 


(2.11) 


aB  =  theoretical  value  of  the  RCS  of  the  reference 
Pt  =  average  power  of  the  target  returns 
PB  =  average  power  of  the  reference  returns 


In  the  following  we  will,  under  certain  assumptions,  attempt  to  estimate 

E 

the  difference  between  the  estimated  target  cross  section  0  and  the 


cataloged  target  cross  section  a  ,  and  determine  statistically  the 


signal -to-noi se  power  ratio. 


We  first  assume  that  the  total  observation  time  does  not  exceed  the 
interval  over  which  the  reference  is  statistically  stationary.  Barrick 
[20]  indicated  that,  for  sea  scatter,  this  total  observation  time  may 
vary  between  1  and  12  hours. 

c 

Suppose  the  target  cross  section  in  the  catalog  can  be  expressed 
as 

oC  =  °B  <pt>  (2.12) 

*  <V 

where  <•>  means  expected  value.  For  convenience  assume  that  the  actual 

value  of  a  can  be  estimated,  and  that  the  target  return  p^  is 
B 

stationary  so  that 

Pt  a  <Pt>  (2.13) 

From  Eqs.  (2.11),  (2.12),  and  (2.13)  we  have 
E 

a  „  aBPt 
t  - 

»  aB<Pt> 

PB 

_  qB<Pt>/<pB> 

V<fV 
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or 


(2.14) 


PB/<PB> 


In  general  we  can  assume  that  the  radar  returns  are  narrow-band  signals 
which  can  be  expressed  as  follows. 


f(t)  =  A(t)cos(u>ct  +  <f(t)) 


(2.15) 


where  is  the  carrier  frequency.  For  convenience  define  a  phasor  f(t) 


as 


f(t)  e  A(t)  eJX^^ 


(2.16) 


Evidently 


*  j  w^t 
f(t)  =  Retf(t)e  } 


(2.17) 


Expand  f(t)  in  terms  of  a  Fourier  series  in  the  interval  [-T/2,  T/2]  and 
obtai n 


j(ka)Ct  +  4>k) 

f(t)  =  Z  Ck  e 
k=-°° 


(2.18) 


where  oi0  =  2tt/T  =  fundamental  frequency,  and  the  Fourier  coefficient  is 
given  by 


j  4>k  T/2  -  -jkb>0t 

Cue  =  1  /  f(t)  e  dt 

T  -T/2 


with  the  variable  being  real.  Substitute  Eq.  (2.18)  into  Eq.  (2.17) 
to  get 


-  6  - 


n 
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f{t)  =  Re{f(t)e  } 


j  ( ki»)0t+<t>|< )  juct 

=  Re(z  Cke  e  } 

k 

-  Z  CkCOS[  (  ojc  +  ka)0)t  +  <j>k] 
k 

=  z  Ck[cos(iuc  +  kco0)t  cos<t>k  _  sin(ojc  +  kw0)t  sin<|>k] 
k 

=  z  [Xk  cos(o>c  +  koi0) t  -  yk  sinfoic  +  ka^H]  (2.19) 

k 

where 


Xk  =  Ck  cos<j>k 


Yk  =  Ck  si n<hk 


The  spectral  component  of  f(t)  at  u»  =  wc  +  ku>0  Is  given  by 


pk  =  A2  +  Yk2  (2.20) 


Assume  that  Xk  and  Yk  are  uncorrelated  Gaussian  variables  with  zero  mean 
2 

and  variance  a  .  This  assumption  is  believed  to  be  applicable  in  most 
k 

cases.  In  fact  Barrick  [20]  experimentally  showed  that  the  spectral 
components  at  the  first-order  sea  scatter  frequency  can  be  approximated 
by  Eq.  (2.20),  with  Xk  and  Yk  being  uncorrelated  zero-mean  Gaussian 
random  variables.  Suppose  M  independent  spectral  samples  of  the  radar 
returns  are  collected,  and  define 


(2.22) 


then  Eq.  (2.14)  can  be  rewritten  as 
c 

E  °t 
a  =  L 

t  T 

Now 

<PB>  =  <XB2  +  YB2> 


=  <XB2>  +  <YB2> 
2 

=  2a 
B 


Let 

M 

PB  =  1  1  PB,n 

W  n=l 

M  2  2 

=  1  z  (*B,n  +  YB,n) 
W  n=l 


(2.23) 


then  r  can  be  rewritten  as 


r 


_PB_ 

<PB> 


^T 

M 

1  Z 
JT  n=l 


M 

Z 

n=l 


B,n 


2 

(YB,n) 

“B 


(2.24 


Since  Xg  n/og  and  Yg  n/aB  are  Gaussian-distributed  with  zero  mean  and 
variance  1,  it  can  be  shown  [31]  that  r  is  Chi -Squared  distributed  with 


M-l  degrees  of  freedom.  The  analysis  will  be  greatly  simplified  if  M  is 
sufficiently  large  so  that  r  can  be  approximately  described  by  a 
Gaussian  distribution  by  virtue  of  the  Central-Limit  Theorem.  It  has 
been  shown  that  M  =  12  would  be  sufficient  [32].  Now  let's  calculate 
the  expected  value  of  r,  namely  <r>,  as  follows. 


<r> 


=  <Pb> 
<Pb> 


Now 


M 

<pB>  =  <  1  E  pB,n  > 
ff  n=l 


.  M  2  2 

=  1  <  Z  X  +  Y  > 

ff  n-1  B»n  B»n 


M 

=  1  z 
TT  n=l 


<X8,n>  +  <YM> 


2 

2a 

B 


where 


2  2  2 

<X  >  =  <Y  >  =  a 
B,n  B,n  B 


(2.25) 


(2.26) 


From  Eq.  (2.23)  and  Eq.  (2.26)  we  have 


<r> 


=  <h> 
<Pb> 
2 
2a 
B 


2ctB 


=  1 


(2.27) 


£ 


The  variance  of  r  can  be  calculated  as  follows. 


Var(r)  =  Var  (  ) 

<Pb> 


=  Var(PB) 


Now 


<v 


Var(PB)  -  Var  (1  r  PB,„I 
M  n=l  ' 


M  2  2 

=  Var  (1  I  X  +  Y  ) 

W  n=l  8,n  B,n 

0  M  2  2 

=  (l)2  Z  Var(X  )  +  Var(Y  ) 
W  n=l  B,n  B,n 


where 


and 


2  4  2 

Var(X  )  =  <XB  n>  -  >2 

B,n  ’  B,n 


4  2  , 

Var  Ye  n  2  =  <Y  >  -  <Y  >2 
°’n  B,n  B,n 
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Therefore  the  variance  of  r  is 

Var(r)  -  <arlpB> 

<PB>2 

4/ 


4  a 
B 


1 

(2.28 


We  conclude  that,  for  M  sufficiently  large,  r  is  approximately 
Gaussian-distributed  with  mean  1  and  variance  1/M.  Let's  rewrite  r  as 

r  =  r0  +  1  (2.29 

where  r0  is  approximately  Gaussian  with  zero  mean  and  variance  1/M.  If 
the  absolute  value  of  r0  is  less  than  1  we  can  rewrite  Eq.  (2.22)  in 
terms  of  a  power  series  expansion  as 


c 


=  <4  (1  -  r0  +  r02  -  r03  +  r04  -  .  .  .  )  (2.30! 

It  is  seen  that  the  estimated  target  cross  section  can  now  be  expressed 

c 

as  the  sum  of  the  cataloged  velue  a  and  a  noise  term  defined  as 


It  can  also  be  shown  that  for  M  greater  than  3  the  probability  that  the 
absolute  value  of  r0  is  greater  than  one  is  extremely  small.  Indeed, 
for  M  =  3, 

Var(r0)  =  }_  =  1/3 
M 

The  probability  that  1 r0 ]  is  smaller  than  1  is  given  as 
p(|r0|  <  1)  =  p{ |3r0!  <  3) 

where  the  term  (3r0)  is  a  standard  Gaussian  variable.  From  the 
Standard  Gaussian  Distribution  table  we  have 

p( |r0|  <  1)  »  99.74%  (2.32) 

The  variance  of  the  noise  N0  can  be  calculated  as  follows. 

Var(N0)  =  at2Var(r0)  +  at2Var(r02)  +  ot2Var(r03)  +  .  .  .  (2.33) 

Eq.  (2.28)  already  gives  us  the  variance  of  r0.  Now  let's  use  Gaussian 
moment  expansion  [21]  to  evaluate  the  higher  order  terms.  From  [21]  we 
have 


„  1  .  3 

E{r0"}  =  { 

0 


(n-1 


’( 


1  n/2 

R  ■' 


for  n  even 
for  n  odd 


Now 


(2.34) 


Var(r02)  =  E { r0^ }  -  E2(r02) 
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(2.36) 


Var(r03)  =  E{rQ6}  -  E2{rQ3}  =  15,u3 

\M/ 

It  is  obvious  that  the  third  and  higher  order  terms  are  negligibly  small 
for  large  M.  Therefore  Eg.  (2.33)  becomes 


Var(N0)  -  (o  )  ( 1  +  ?  .  ,'1\  ) 

M  \MV 


(2.37! 


Define  a  signal -to-noise  power  ratio  as 


c 

a 


(2.38] 


From  Eg.  (2.37)  and  Eg.  (2.38)  we  have 


As  an  example  let's  suppose  the  coherent  integration  time  for  each 
sample  is  10  seconds,  and  we  take  M  =  100  samples,  which  corresponds  to 
a  total  observation  time  of  about  16  minutes,  then  the  post-processing 
signal-to-noise  ratio  is 


s  ,  100 
tt  "ynr? 


10  =  10  dB 


In  Chapter  V  we  will  determine  the  probability  of  mi sclassifi cation 
corresponding  to  this  observation  period. 


D.  MEASUREMENT  OF  PHASE  RETURNS 

In  this  section  we  shall  briefly  analyze  the  phase  returns  measured 
by  a  typical  coherent  radar.  Fig.  2.7  is  a  simplified  block  diagram  of 
a  coherent  pulse  radar  system.  Suppose  the  RF  oscillator  generates  a 
sinusoidal  waveform  of  frequency  fRp,  and  assume  that  the  output  signal 
is  given  by 

E  =  I E0 1  cos(2irfRpt)  (2.40) 

The  signal  transmitted  by  the  antenna  can  be  given  as 

Ey  =  |Ey|  cos(2TrfRpt  +  <|>y)  (2.41) 

where  <j>y  is  the  total  phase  shift  between  the  RF  oscillator  and  the 
antenna.  The  received  signal  at  the  RF  amplifier  of  the  receiver  will 
then  be 


39 


Er  =  | Er I  cos ( 2iTf Rpt  +  <t>  +  <j>q  +  4-ttR  )  (2.42) 

where  4>r  is  the  total  equipment  phase  shift  between  the  RF  oscillator 
and  the  input  to  the  receiver  measured  at  the  equipment  operating 
frequency,  4>  is  the  intrinsic  phase  return  induced  by  the  target,  and  R 
is  the  one-way  wave  path  length.  The  local  oscillator  frequency  is 
obtained  by  mixing  the  output  of  the  RF  oscillator  with  the  output  of  a 
stable  IF  reference  oscillator  and  then  selecting  either  the  sum  or 
difference  frequency  with  a  bandpass  filter.  The  local  oscillator 
signal  is  then  mixed  with  the  output  of  the  RF  amplifier  to  obtain  the 
IF  frequency,  which  is  in  turn  amplified  by  the  IF  amplifier.  At  this 
point  the  signal  is  of  the  form 

Eif  =  I E IP |  cos ( 2irf ipt  +  41  +  4ttR  +  4>jp)  (2.43 

X 

where  d>ip  is  the  total  equipment  phase  shift  at  the  output  of  the  IF 
amplifier.  The  output  of  the  amplifier  is  coherently  demodulated,  with 
both  in-phase  and  quadrature  signals  from  the  IF  reference  oscillator 
used  as  a  reference. 

This  results  in  the  following  in-phase  and  quadrature  components: 

EI  =  1  ElF l  cqs(47tR  +  4>  +Hotal )  (2.44 

x 

eQ  =  !  E1F 1  sin(4wR  +  +  Hotal) 


where  4>total  1S  the  total  system  phase  shift.  The  phase  output  of  the 


signal  processor  is  given  by 


0  =  tan-1  (£0.) 

^  I 

=  4^tR  +  4>  +  4>total  (2.45) 

X 

The  wave  path  length  R  for  sky  wave  radars  is  a  function  of  ionospheric 
condition  and  frequency,  and  the  phase  may  vary  widely  as  the  frequency 
changes.  For  surface  wave  radars  the  wave  path  length  is  essentially 
independent  of  frequency. 

Goggins  and  others  [19]  showed  that  it  would  be  possible  to  recover 
the  intrinsic  phase  returns  under  certain  conditions.  They  [19] 
mathematically  proved  that  for  a  smooth  conducting  body  in  the  Rayleigh 
region  the  co-pol ari zed  terms  of  the  scattering  matrix  have  zero  phase 
shift.  Based  on  this,  they  suggested  a  differential  phase  technique  to 
recover  the  intrinsic  phase,  as  discussed  below. 

From  Eq.  (2.45)  the  phase  return  measured  in  the  Rayleigh  region 
can  be  given  as 

60  =  4wR  +  <j>°  (2.46) 

total 

where  x0  is  the  wavelength  in  the  Rayleigh  region,  R0  is  the  wave  path 
length  in  the  Rayleigh  region,  and  the  intrinsic  phase  <|>  is  set  to  zero. 
The  total  equipment  phase  shift  can  be  written  as 

♦total  =  (2-47) 

X 


where  D  =  Total  wave  path  length  In  the  equipment.  From  Eq.  (2.45)  and 
Eq.  (2.46)  we  have 

e  =  4ttR  +  ,(,  +  2*D 

Qo  =  47r^o  + 

^o  ^o 
or 

X0  =  4ttR  +  X  <|>  +  2nD  (2.48) 

XO0O  =  4tvR0  +  2hO  (2.49) 

Taking  the  difference  between  Eq.  (2.48)  and  Eq.  (2.49),  we  have 
X0  -  XO0O  =  4tt(R-R0)  +  X<j> 

If  R  =  R0  then  the  intrinsic  phase  is  given  by 

$  =  A9  ~  xo9o  (2.50) 

X 

For  sky-wave  radars  the  wave  path  length  R  is  obviously  not  equal  to  R0, 
and  Eq.  (2.50)  does  not  hold.  For  surface-wave  radars,  R  may  be  close 
to  R0,  and  the  intrinsic  phase  may  be  recovered.  In  addition  to  the 
requirement  that  R  be  equal  to  Rq,  this  method  needs  phase  measurements 
in  the  lower  Rayleigh  region,  which,  for  targets  like  ships  and 
aircraft,  corresponds  to  frequencies  much  lower  than  the  H.F.  region. 
Because  of  such  requirements,  the  classification  alogorithms  to  be 
discussed  in  the  next  chapter  do  not  directly  utilize  the  intrinsic 
phase. 
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CHAPTER  III 


TARGET  CLASSIFICATION  ALGORITHMS 

A.  INTRODUCTION 

The  subject  of  pattern  classification  has  been  extensively  studied 
in  the  literature  (e.g.,  Duda  and  Hart  [21]).  The  basic  process  of  a 
pattern  classification  system  is  depicted  in  Figure  3.1.  As  shown  in 
the  figure,  the  measurement  system  (in  our  case  a  radar  system)  makes 
measurement  of  a  physical  object,  and  each  measurement  may  contain  a  set 
of  independent  components  and  can  be  viewed  as  a  vector  y.  In  our  case 
each  measurement  is  a  vector  whose  components  are  the  amplitude  and 
phase  of  the  radar  returns  of  a  target  at  a  number  of  frequencies.  Each 
measurement  vector  is  a  point  in  the  multi -dimensional  measurement  space 
The  purpose  of  the  feature  extractor  is  to  reduce  the  dimensionality  of 
the  measurement  vector  by  using  certain  features  that  distinguish  the 
objects.  The  feature  extractor  transforms  a  measurement  vector  y  to  a 
point  or  feature  vector  in  a  mul ti-dimensional  feature  space.  The 
feature  vectors  are  passed  to  the  object  classifier  that  uses  the 
features  and  assigns  the  object  to  an  appropriate  class. 


MEASUREMENT  FEATURE 

VECTOR  VECTOR 


Figure  3.1.  Basic  process  of  a  pattern  classification  system. 
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There  are  two  types  of  classification  problems,  i.e.,  parametric 
and  nonparametric.  In  parametric  classification  the  probability 
distribution  of  the  observation  space  is  expressed  as  a  functional  form, 
such  as  that  of  the  normal  distribution,  which  is  completely  specified 
by  at  most  a  finite  number  of  real -valued  parameters  [33].  In 
nonparametric  classification,  on  the  other  hand,  the  probability 
distribution  is  not  indexed  by  a  finite  number  of  parameters  [33].  Note 
that  the  fundamental  difference  of  the  nonparametric  classifier  from  the 
parametric  classifier  is  that  the  parametric  form  of  the  underlying 
probability  distribution  is  not  assumed  or  used.  The  nonparametric 
method  are  important  because  parametric  probabilistic  information  is 
usually  not  available  or  the  assumption  of  it  is  not  justified. 

In  Section  B  we  will  briefly  discuss  a  fundamental  classification 
theory,  i.e.,  the  Bayes  decision  theory  [21,22].  As  will  be  seen,  such 
decision  theory  minimizes  the  average  probability  of  error,  but  requires 
knowledge  about  the  a  posteriori  probability  of  the  feature  vector  X. 

For  classification  of  radar  targets,  it  is  clear  that  parametric 
probabilistic  information  about  the  feature  vector  is  not  available  or 
the  assumption  of  it  is  not  justified  because  of  the  lack  of  adequate 
information  on  system  performance.  Therefore,  we  must  resort  to 
nonparametric  methods.  Two  classification  methods  which  do  not  require 
knowledge  of  the  probabilistic  information  will  be  discussed  in  Section 
C  (the  frequency-domain  algorithm)  and  Section  D  (the  time-domain 
algorithm). 
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8.  BAYES  DECISION  THEORY 

Suppose  two  classes  of  samples  are  to  be  classified.  Let  the  first 
class  be  denoted  as  and  the  second  as  Let  the  probability  of 
occurrence  (a  priori  probability)  of  class  w-j ,  i=1.2,  be  denoted  as 
p(u>i).  Let  p{ X/ui-j )  denote  the  conditional  probability  of  an  observation 
X  in  the  feature  space,  given  that  X  belongs  to  class  w-j .  Given  an 
observation  X,  the  probability  that  the  observation  belongs  to  class  w-j 
is  given  as  piu^/X),  i.e.,  the  a  posteriori  probability.  Suppose  that 
we  know  both  the  a  priori  probabilities  p(o>j)  and  the  conditional 
densities  p(wj/X),  then  the  a  posteriori  probability  can  be  calculated 
from  the  Bayes  rule: 

p( X/ui-j)  p(u)-j) 

PH/X)  =  -  (3.1) 

p(x) 

where 

2 

p(  x )  =  l  p(x/o)j)  p(wj)  . 

3  =  1 

Because  we  have  assumed  that  the  decision  rule  must  classify  an 
observation  X  either  as  class  1  or  class  2,  we  can  view  it  as  a  rule  for 
dividing  the  total  observation  space  Z  into  two  parts,  Z^  and  Z2.  When¬ 
ever  an  observation  falls  in  1\  it  is  classified  as  class  1,  and  when¬ 
ever  an  observation  falls  in  Z2,  it  is  classified  as  class  2.  We  can 
now  write  an  expression  for  the  average  cost  or  risk  in  terms  of  the 
decision  regions 
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RISK  =  P(Wi )  /  P(X/Wi)dX  +  P(W2)  /  P(X/W2)dx  (3.2) 

Z2  Zi 

where  we  have  assumed  that  there  is  no  cost  for  making  right  decisions, 
and  the  costs  for  making  wrong  decisions  are  unity.  Now,  to  find  the 
Bayes  decision  rule  we  must  choose  the  decision  regions  1\  and  Z2  in 
such  a  manner  that  the  risk  will  be  minimized  [21,22].  Since 

Z  =  Zi  +  Z2 

Eq.  (3.2)  can  be  rewritten  as 

RISK  =  P(W!)  /  P(X/Wi)dx  +  P(W2)  /  P(X/W2)dx 
Z-Zl  Zi 

or 

RISK  =  1  +  /  [P(W2)  P(X/W2)  -  P(Wi)  P(X/Wi )]  dx 

Zl 

In  order  to  minimize  the  risk  all  values  of  X  where  the  second  term  in 
the  above  integral  is  larger  than  the  first  should  be  included  in  1\ 
because  they  contribute  a  negative  amount  to  the  integral.  Similarly, 
all  values  of  X  where  the  first  term  is  larger  than  the  second  should  be 
excluded  from  Z^  (assigned  to  Z2)  because  they  would  contribute  a  posi¬ 
tive  amount  to  the  integral.  Values  of  X  where  the  two  terms  are  equal 
have  no  effect  and  may  be  assigned  arbitrarily.  Thus,  the  decision 
regions  are  defined  by  the  statement: 

If  P(W2)  P(X/W2)  >  P (Wi )  P(X/Wi) 
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assign  X  to  Z2  and  consequently  classify  X  as  class  2.  Otherwise  assign 
X  to  1\  and  classify  X  as  class  1.  In  terms  of  the  a  posteriori  prob¬ 
ability,  the  Bayes  decision  rule  is  given  as  follows. 

classify  X  as  u>i  if  p{a^/X)  >  pfa^/X) 

classify  X  as  w 2  if  p(u>2/X)  >  p( uj^/X )  .  (3.3) 

If  there  are  more  than  two  classes,  say,  N  classes,  the  Bayes  decision 
rule  becomes  [21,22] 

classify  X  as  if 

p ( Lo-f /X )  >  p(uj/X)  for  j  =  1,2,...,  N  and  i  *  j  .  (3.4) 

It  is  evident  that  the  Bayes  decision  rule  requires  knowledge  about  the 
a  posteriori  probability  density  function  which  is  not  available  for 
application  to  this  problem.  In  the  following  two  sections  two  classi¬ 
fication  algorithms  which  do  not  require  the  probabilistic  information 
will  be  discussed. 

C.  CLASSIFICATION  IN  THE  FREQUENCY  DOMAIN 

The  first  classification  algorithm  makes  use  of  features  extracted 
from  the  multi -frequency  radar  returns  of  a  target.  The  measurement 


vector  may  be  chosen  as  a  set  of  amplitude  A(fj)  and  phase  e ( f ■} )  of  the 


backscatter  of  a  target,  measured  at  frequencies  f  =  fj.,  f2,  ....  fuj. 
Note  that  the  amplitude  is  defined  as  the  square  root  of  the  measured 
target  cross  section.  From  Equation  (2.45)  we  know  that  the  relative 
phase  return  0  is  a  function  of  the  wave  path  length  R  and  the  total 
equipment  phase  shift  <l>total  ♦  i*e., 


9  =  4ttR  +  $  +  <j>.  .  , 

total 


(3.5) 


As  discussed  in  Chapter  II,  section  D,  the  intrinsic  phase  return  can 

be  recovered  by  the  differential  phase  technique  when  the  phase  path 

length  R  is  independent  of  the  frequency.  The  technique  also  requires 
that  a  phase  return  be  measured  in  the  Rayleigh  region. 

Instead  of  using  either  the  relative  phase  return  0  or  the 

intrinsic  phase  return  <t>,  we  can  define  a  "differential  quantity"  W  as 

follows: 


From  Equation  (3.5)  we  have 


0.X.  =  4ttR.  +  <{>.  A.  +  ^0tal  X. 
11  1  1  M  1 


9i+lXi+l  =  4irRi+l  +  *i+lx1+l  +  ‘•’i+l31  Xi+1 


(3.6) 


where  0j  =  measured  phase  return  at  x  =  Xj 

R-j  =  one-way  wave  path  length  at  x  =  Aj 
4>j  =  intrinsic  phase  at  X  =  A,- 


Except  for  the  factor  4ir(Rf-Ri+i) ,  the  quantity  ( e-j )  is 
directly  related  to  the  quantity  UiM'^l+lM+l) ,  which  is  a  function  of 
only  the  wavelengths  and  the  intrinsic  phases.  For  surface  wave  radars, 
we  have  R-j  »  Ri+i,  so  that 

®i  Xj  -  ®i+lXi+i  "  <H  Xj  -  'H+lX-j+1.  (3.9) 

For  sky  wave  radars,  R-j  and  R-f may  differ  by  a  large  amount.  However 
if  x,  and  Aj+i  differ  by  only  a  small  amount  and  the  radar  operating 
frequencies  are  properly  chosen,  the  difference  between  R-j  and  Rj+i  may 
be  made  very  small.  Define  W  as 

Wi  =  9j  A-j  -  9^+i  A1+i 

1  *  1,  2,  ....  N-l  (3.10) 


m 


■ 


51 


where  N  is  the  total  number  of  frequency  components  available.  Since 
the  phase  0j  in  Equation  (3.5)  is  multiply  valued,  i.e., 

0j  =  0-j  +  2htt  ;  n  *  1,  2,  3,  ...,  (3.11) 

the  variable  W-j  in  the  equation  has  also  multiple  values.  This  problem 
can  be  solved  as  follows. 

Let  the  phase  difference  A0  be  defined  as 

A0  =  0i  -  01+1  (3.12) 

It  is  evident  that  the  variable  W,  is  uniquely  defined  if  we  require 
that  the  phase  difference  be  confined  to  the  region: 

<  A0  <  TT  (3.13) 

Since 

A®  =  9j  -  Qi+1 

=  (  Qi  -  9i+i)  +  2n7r  ;  n  =  integer 
the  phase  difference  A0  can  be  confined  to  the  region  defined  in 
Equation  (3.13)  by  repeatedly  subtracting  2v  or  adding  2tt  to  A9  until 
Equation  (3.13)  is  satisfied. 

Define  a  feature  space  as 

(Ai,  A£,  ...,  An,  kWi,  kW2,  ...»  kW|q_i )  (3.14) 

where  K  is  a  normalization  constant  such  that  the  variances  of  A  and 
W  are  related  as 


VAR(A)  =  0  •  VAR(kW) 


(3.15) 


If  B  Is  set  to  1  then  both  the  amplitude  and  the  phase  are  assigned 
equal  weightings.  If  the  phase  is  not  reliable  or  contains  large  amounts 
of  error,  as  is  the  case  for  ionospheric  propagation,  we  can  set  3  to  be 
a  large  number  so  that  the  amplitude  is  given  more  weight.  On  the  other 
hand  3  can  be  set  to  be  a  small  number  so  that  the  phase  is  given  higher 
weighting.  Thus,  by  setting  3,  we  can  determine  the  reliable  weights  of 
the  amplitude  and  phase.  Note  that  if  3  is  set  to  infinity,  then,  only 
the  amplitude  is  used  for  classification,  and  if  3=0,  then,  only  the 
phase  is  used  for  classification.  Once  3  is  selected  the  normalization 
constant  K  can  be  computed  as 

(3.16) 

The  problem  now  is  how  to  estimate  the  variance  of  the  amplitude  and  W. 
Our  approach  here  is  to  calculate  the  variance  of  A  and  of  W  in  the 
data  base,  and  use  the  result  as  an  estimate  of  the  true  variance.  It 
is  noted  that  such  estimates  may  contain  errors,  and  further  research 
into  finding  an  optimum  K  is  needed. 

Once  the  constant  K  is  calculated,  the  distance  between  any  two 
points  in  the  feature  space  Equation  (3.14)  can  be  defined  as 


where 

Aa  and  Wa  define  one  of  the  two  points  and,  Ab  and  Wb  define  the  other 
point. 
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^ ■  V’  v  v  t 


"T*T 1 


Once  the  feature  space  is  defined,  the  Nearest  Neighbor  Decision 
rule  [21]  can  be  applied  to  classify  the  targets.  The  N.N.  rule  is  non- 
parametric,  yet  its  error  rate  is  close  to  the  error  rate  of  the  Bayes 
decision  rule.  In  fact,  it  has  been  shown  [21]  that  the  error  rate  of 
the  N.N.  rule  is  never  worse  than  twice  the  Bayes  error  rate.  An  algor¬ 
ithm  for  classification  in  the  frequency-domain  is  given  below. 

Suppose  we  choose  the  frequencies  f  =  f^,  f?,  ...,  f^.  Let  the 
amplitude  and  phase  returns  of  the  unknown  target  be  given  as  A^( fi )  and 
0t(f-j),  respectively,  where  i  =  1,  2,  ...,  N.  Also  suppose  there  are  M 
classes  in  the  data  base.  Let  the  amplitude  and  phase  returns  of  each 
class  in  the  database  be  given  as  Aj(f-j)  and  8j ( f ) ,  where  fi  means  the 
ith  frequency,  i  =  1,  2,  ...,  N,  and  j  means  the  jth  target,  j  =  1,  2, 

. . .  ,  M . 


1. 


Compute  the  W's  as  follows 

“f  ■  Vt<v  -  ViVfm 
=  »i VV  - 


) 

) 


i  =1,  2,  ...,  N-l 
j  =  1,  2,  ... ,  M 


2.  Calculate  the  means  of  A  and  W  in  the  database. 
M  N 

J  =  1  z  l  Aj  (f,) 

WT  j=l  i=l 


M  N-l  j 

¥  =  1  E  E  Wi 

MTTT-TT  j  =  l  i  =  l 
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D.  CLASSIFICATION  IN  THE  TIME  DOMAIN 


In  the  last  section  we  tried  to  reduce  the  phase  error  caused  by 
the  variation  of  the  effective  phase  path  length.  Since  this  phase 
error  is  equivalent  to  time-shifting  of  the  corresponding  time-domain 
waveform,  it  may  also  be  removed  in  the  time  domain.  It  was  shown  [12] 
that  the  time-domain  ramp  response,  obtained  by  integrating  twice  the 
impulse  response,  is  proportional  to  the  cross-sectional  area  of  a 
target  along  the  line-of-sight.  Here  the  impulse  response  waveform  is 
the  time-dependent  electro-magnetic  field  intensity  produced  at  a  fixed 
point  in  space  when  a  plane  electromagnetic  shock  wave  strikes  an  object 
[3],  and  can  be  obtained  by  synthesizing  the  frequency-domain  response 
[12].  It  is  therefore  feasible  to  use  the  time-domain  data  as  features 
for  identification  of  targets. 

Let  X(m)  be  a  set  of  complex  frequency-domain  data  defined  by 
j  e(fk) 

X(m)  =  A(fk)  e  *,  k  =  1,  2,  ...,  M  (3.20) 

where  A  is  the  amplitude  and  0  is  the  phase  of  the  target  returns 

m  =  fl/ Af  >  VAf  >  •••»  WAf  (3.21) 

and  the  frequencies  are  equally-spaced,  i.e., 

fi  '  fi-l  =  Af  ;  1  =  1,  2,  ...,  M  (3.22) 

Define  the  time-domain  data  as  the  Discrete  Inverse  Fourier  Transform 
(DIFT)  of  the  frequency  domain  data  X(m)  [34],  i.e., 
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(3.23) 


(3.24) 


Note  here  the  time-domain  data  should  be  real,  and  the  amplitude  of  the 
frequency-domain  data  are  symmetric  around  f  =  0,  i.e., 

X (m)  =  X*(-m)  (3.25) 


where  X*  is  the  complex  conjugate  of  X.  Suppose  the  time-domain  data  of 
two  targets,  x(n)  and  y(n),  respectively,  are  synthesized  from  their 
respective  frequency-domain  response.  Define  a  Correlation  Function  of 

x(n)  and  y(n)  as 
N 

R  Jk)  =  E  x(n)  y(n-k)  (3.26) 

xy  n=l 


From  the  Cauchy-Schwartz  inequality,  i.e., 

I  Ea-jbj  1 2  <  (E|ai|2)  (E|bi  |2)  (3.27) 

i  i  i 

we  have 

N  N  N 

[Rxy(k)]2  =  [  Z  x(n)  y ( n-k ) ] 2  <;  (  z  x2(n))  (  z  y2(n))  (3.28) 

n=l  n=l  n=l 
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where  the  equality  holds  if  and  only  if 


Define  a  correlation  coefficient  of  x{n)  and  y(n)  as 


pXy(k)  = 


Rxy(k) 

Usn^{ 


It  is  evident  from  Equations  (3.28)  and  (3.32)  that 


(3.33) 


Pxy ( k )  <1 


(3.34) 


and  the  equality  holds  if 


x=y  and  k=0 


(3.35) 


i.e.,  if  the  two  waveforms  coincide  and  the  relative  time-shift  between 
the  two  waveforms  is  zero. 

A  time-domain  classification  algorithm  is  given  below.  Suppose  the 
time-domain  data  x(n)  of  a  target  is  derived  from  its  frequency-domain 
response,  and  the  time-domain  data  of  each  class  in  the  data  base  are 
derived  from  their  respective  frequency - doma i n  data,  and  are  given  as 
y-f  ( n)  ;  1=1,2,...,  C. 


1.  Compute  the  Correlation  Function 
R  (k) 

Pxy.(k)  x*i _ 

1  /Wnri 


1  =  1,  2,  ....  C 


(3.36) 
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2.  Choose  the  time-shift  constant  K  such  that 
Pxy^ (k)  is  maximized,  i.e.. 


pxy-  =  MfX  t  pxyi(k)}  for  each  1 


(3. 


Note  that  this  would  remove  the  relative  time-shift  between  the  two 

waveforms,  and  p™a*  is  independent  of  the  time-shift, 

xyi 

3.  Identify  the  unknown  target  as  class  i  if 


max 

D 

xyi 


>  p 


max 

xyj 


j  =  1,  2, 


j  *  i 


(3. 


Note  that  from  the  theories  of  Discrete  Fourier  Transform  [34] 


Rxy(k)  =  DIFT  {x(m)  y*(m)} 


(3. 


and 


N  M 

Rxx(0)  =  z  X2(n)  =  2  z  lX(m)|2  (3. 

n=l  m=l 


Thus  it  is  more  efficient  to  calculate  the  correlation  coefficient  as 


CHAPTER  IV 


GENERATION  OF  A  DATA  BASE 

A.  INTRODUCTION 

The  testings  of  classification  algorithms  require  a  data  base  which 
is  representative  of  the  data  to  be  expected  in  practice.  A  set  of 
radar  backscatter  measurements  made  by  The  Ohio  State  University 
ElectroScience  Laboratory  of  scaled  models  of  ships  and  aircraft  was 
used  to  generate  a  data  base  for  testing.  The  data  base  was  designed  to 
simulate  a  full  scale  data  base  for  target  classification  using  H.F . 
multi  frequency  radars. 

In  section  B  we  will  discuss  the  measurement  of  scaled  models  at 
scaled  frequency  and  the  calibration  of  measured  radar  returns.  Due  to 
limitations  on  the  measurement  system,  backscatter  of  the  models  was 
measured  in  a  number  of  small  frequency  bands  respectively,  e.g.,  2-4 
GHz,  4-8  GHz,  and  8-12  GHz.  Target  locations  may  be  slightly  changed 
from  one  measurement  to  another,  therefore,  when  these  data  were 
connected  together  to  form  a  wider  band  of  data,  e.g.,  2-12  GHz,  phase 
offsets  were  found  to  be  present  at  the  junction  of  two  adjacent 
frequency  bands.  In  section  C  we  will  discuss  the  connection  of  two 
data  sets  and  techniques  for  removing  the  phase  offsets.  After 


cal i oration  and  phase  correction,  the  data  were  interpolated  and  smooth¬ 
ed  to  obtain  data  values  at  desired  frequencies.  The  interpolation 
and  smoothing  made  use  of  a  Hamming  window  which  had  zero  phase  in  both 
the  frequency  domain  and  the  time  domain  [35],  and  was  equivalent  to  a 
low-pass  filter  in  the  time  domain.  The  interpolation  will  be  the  sub¬ 
ject  of  Section  D. 

In  summary,  the  measured  data  were  first  calibrated,  then  connected 
together  in  bands  and,  then,  interpolated  and  smoothed.  The  resulting 
data  formed  the  data  base  which  was  used  for  testing  the  algorithms 
discussed  in  Chapter  III.  In  section  E  we  will  present  plots  of  the 
post-processed  waveforms  in  the  data  base. 

B.  MEASUREMENT  AND  CALIBRATION 

The  measurement  facility  used  to  obtain  scattering  data  on  the 
model  targets  was  a  monostatic  swept-frequency  radar  system.  Operation 
of  the  system  was  interactively  controlled  by  a  FORTRAN  program  on  a 
minicomputer.  A  simplified  diagram  of  the  system  is  shown  in  Fig.  4.1. 
Model  targets  were  placed  on  a  low-RCS  pedestal.  The  ranqe  was 
configured  so  that  the  far-field  criterion  was  approximately  satisfied. 
The  elevation  angles  of  the  incident  waves  were  approximately  zero.  For 
ship  models  a  large  ground  plane  is  needed  to  simulate  the  reflection 
from  the  sea  surface.  Instead  of  buildinq  such  a  ground  plane,  a 
replica  was  attached  to  each  ship  model,  as  shown  in  Fiq.  4.2.  This 
effectively  simulated  the  image  created  by  the  reflection  from  the  sea 
surface.  By  rotating  the  targets,  the  aspect  angles  of  the  incident 


A  SHIP  MODEL 


A  REPLICA 
SHIP  MODEL 


OF  THE 


Fig.  4.2.  Simulation  of  the  reflection  from  the  sea  surface, 
replica  is  attached  to  each  ship  model. 


waves  were  varied.  For  each  particular  aspect  angle,  the  radar  returns 
were  measured  in  a  number  of  frequency  bands,  which  included  1-2  GHZ, 
2-4GHZ,  4-8GHZ,  and  8-12GHZ.  In  each  band  there  were  201  equally-spaced 
points.  The  frequency  increments  of  the  data  points  were  0.005GHZ  for 
the  1-2GHZ  band,  0.01GHZ  for  the  2-4GHZ  band,  .02GHZ  for  the  4-8GHZ 
band,  and  .02GHZ  for  the  8-12GHZ  band.  The  frequency  increments  were 
chosen  to  be  small  to  provide  enough  data  points  for  the  interpolation 
and  smoothing  of  the  data  waveform. 

The  aspect  angles  and  the  frequency  bands  chosen  for  each  model 
might  be  different.  In  most  cases,  the  aspect  angles  varied  from  0  to 
180  degrees  at  15  degree  increments,  and  the  frequency  bands  were  chosen 
as  2-4GHZ,  4-8GHZ,  and  8-12GHZ.  In  most  cases,  vertical  polarization 
was  used. 

The  calibration  of  measured  data  made  use  of  reference  spheres 
whose  theoretical  radar  returns  were  readily  available.  Prior  to  cali¬ 
bration  by  spheres,  background  returns  were  subtracted  from  the  data. 

The  background  returns  were  measured  separately  after  the  targets  were 
removed.  They  included  returns  from  the  antenna's  supporting 
structures,  the  building's  walls  and  ceilings,  wave  absorbers,  etc.  It 
was  noted  that  the  subtraction  of  background  returns  would  not  remove 
couplings  between  targets  and  the  background  itself.  Therefore  a 
Hamming  Window,  which  is  equivalent  to  a  time-domain  low-pass  filter, 
was  applied  to  alleviate  the  coupling  effects.  The  voltage  amplitude  of 
the  received  target  return  was  proportional  to  the  square  root  of  the 
target  cross  section.  The  received  voltage  amplitude  A  and  phase  0  for 
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a  particular  frequency  f  =  fj  formed  a  complex  number  given  by 

s(f-j)  =  A( f ^  )cos e{ f i )  +  j  A(f-j)sine(f0) 

The  radar  return  of  a  target  at  f  =  fi  was  calibrated  as 

TCW  =  (Tm(fi)  -  Bm(f1))  •  St(fj) 

^ni(^i)  -  Bm(f-j) 

where  all  items  were  complex  numbers,  and 

Tc(fi)  =  calibrated  radar  return  of  the  target 
Tm(fi)  =  measured  radar  return  of  the  target 

Bm(fj)  =  measured  radar  return  of  the  background 

Sm(f-j)  =  measured  radar  return  of  a  sphere 

St ( f i )  =  theoretical  radar  return  of  the  sphere 

Eq.  (4.2)  is  a  calibration  equation  which  employs  the  backscatter  from  a 
sphere  as  a  reference,  and  is  similar  in  principle  to  the  calibration 
equation  (2.8  ),  which  employs  the  backscatter  from  the  sea  surface  as  a 
reference.  The  phase  reference  of  the  calibrated  return  of  the  target 
is  the  same  as  that  of  the  theoretical  return  of  the  sphere,  i.e.,  the 
center  of  the  target. 

C.  CONNECTION  OF  TWO  FREQUENCY  BANDS 

This  section  is  concerned  with  the  connection  of  two  sets  of  data 
measured  in  two  adjacent  frequency  bands  to  create  a  data  set  of  a  wider 


(4.1) 


(4.2) 


bandwi dth. 


Suppose  two  sets  of  multi -frequency  noise-contaminated  radar 
returns  measured  at  adjacent  frequency  bands  are  to  be  connected  as 
shown  in  Fig.  4.3.  Let  the  frequency  range  of  the  first  set  of  data  be 
f_l  <  f  <  f0  and  that  of  the  second  set  of  data  be  f0  <  f  <  fj.  The 
junction  of  the  two  sets  is  at  f  =  f0.  Suppose  the  positions  of  the 
target  in  these  two  measurements  differ  by  only  a  small  amount.  Note  a 
2-mil imeter  position  offset  at  a  frequency  of  6  GHZ  yields  almost  30 
degrees  of  electrical  phase  discontinuity  (Eq.  (4.4)).  Since  the  target 
is  assumed  to  be  in  the  far-zone,  there  is  no  significant  amplitude 
offset  at  f  =  f0.  Therefore,  the  amplitude  data  can  be  connected 
directly  without  correction.  The  phase  offset,  however,  is  significant 
even  if  the  position  offset  is  only  a  fraction  of  a  wavelength.  The 
position  offset  causes  a  ramp- type  offset  in  phase  as  the  frequency 
increases  or  decreases.  This  type  of  offset  can  be  approximately 
corrected  as  follows. 

Let  the  phase  offset  A*e(f0)  at  f  =  f0  be  estimated  as 

A*e(f0)  =  e  (f0)  -  e*(f0)  (4.3) 

*  ★ 

where  ^(fo)  and  9^^o^  are  estimates  of  the  true  phase  02(fo)  and 
0l(fo),  respectively.  Since  the  phase  offset  and  the  position  offset 
6  are  related  by 

A9(f0)  =  2*-26  ;  X0  =  c/f0  (4.4) 

~ ^o 

The  position  offset  5  can  be  estimated  as 


67 


Now  if  f/f0  is  not  an  integer  02c(f)  has  different  values  for  different 
values  of  n.  A  principal  branch  must  therefore  be  chosen  for  A*e(f0). 

A  convenient  choice  would  be 


-IT  <  A*  0  (  f  Q  )  <  It 


(4.10) 


In  order  to  calculate  A*e(f0)  two  parameters  9*(f0)  and  0* ( f0 )  must  be 

estimated  first.  One  approach  is  to  interpolate  0*(f)  and  0*(f)  using 

1  2 

data  points  in  the  vicinity  of  f  =  f0.  The  minimum-squared-error  linear 
interpolation  can  be  chosen  to  estimate  0*(fo)  and  0*(fo)  as  follows 
[36]. 


e*(f0)  =  Ai f0  +  Bi ;  i  =  1,2 


(4.11) 


where 


e  0i(fj)(fj-n 
iii _ 

"  (fi  -  T)2 


B.  =  9.-  A.? 
1  I  1 


f  -  i  s  f j 

n  j  =  l 


9;=i  *  8i(fj) 
n  j=l 


(4.12) 


-  *  ■  <»  •*  *  "  H  -  w  ■■ ' 

'1  -  s 


L  JW  -*■  L.JL  —  n  — ! 


An  empirical  value  for  n  is  5  or  6.  A  higher-order  interpolation  can  be 
used  so  that  it  fits  better  the  data  points.  However  it  would  likely  be 
distorted  by  the  noise  already  contained  in  the  data. 

The  disadvantage  of  the  interpolation  is  that  a  few  of  the  data 
points  used  might  be  so  enormously  distorted  by  noise  that  the 
interpolation  results  become  unacceptable.  A  remedy  to  this  dilemma  is 
to  estimate  the  phase  offset  A*e(f0)  in  the  time  domain.  A  ramp-type 
phase-shift  in  the  frequency  domain  is  equivalent  to  a  constant 
time-shift  in  the  time  domain.  If 
F 

m(t)  - - -  A(f)e-i0(f>  (4.13) 

form  a  Fourier  transform  pair  then 
F 

m(t-t)  - - -  A( f)ej(  9(f)-2irfT) 

From  Eq.  (4.5)  we  have 

T  -  A*e(f0) 

TttTo” 

or 

A*0(fo)  =  2wf0T  (4.16) 

Suppose  the  time-domain  waveform  is  given  by 
m(t)  =  F-l {A( f )eJ  9< f > } 
and 


(4.14) 


(4.15) 
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(4.17) 


mi(t)  =  F-l{A1(f)e  } 

“  -  j  ©2(  f) 

m2(t)  =  F-l(A2(f)e  } 

where  mi(t)  and  m2(t)  are  band-limited  approximations  to  the  time  domain 
waveform  m(t)  using  different  frequency  bands.  A  plot  of  m^t)  and 
m2(t)  is  shown  in  Fig.  4.4.  The  amount  of  time-shift  t  can  be 
approximated  in  two  ways.  The  first  is  to  roughly  estimate  the  distance 
between  the  peaks  of  mj(t)  and  m2(t).  The  second  is  to  cross-correlate 
the  two  waveforms,  and  find  the  time-shift  which  maximizes  the 
cross-correlation  function.  The  time  approach  uses  the  contribution 
from  all  available  frequency  data  points,  and  the  overall  shapes  of  the 
time  domain  waveforms  cannot  be  adversely  affected  by  a  few  distorted 
data  points  in  the  frequency  domain.  A  drawback  is  the  additional 
processing  time  required.  It  was  decided  that  the  frequency-domain 
linear  interpolation  approach  would  be  appropriate  for  our  application. 
Therefore  the  frequency-domain  approach  was  used. 


D.  INTERPRETATION  AND  SMOOTHING 

In  this  section  we  will  discuss  the  use  of  Hamming  Windows  for 
interpolating  and  smoothing  data  points  to  estimate  the  values  of  the 
data  at  selected  frequencies.  Suppose  the  data  points  are  a  set  of 
complex  numbers  c(f)  given  as 

c(fk)  =  X(fk)  +  j  Y(fk)  (4.18) 


where 


k  =  1,  2,  .  .  .  ,  L 

The  real  part  X  and  the  imaginary  part  Y  are  real  functions  of 
frequency,  as  shown  in  Fig.  4.5.  Define  a  Hamming  Window  as  follows 
[24] 


H( f , t)  =  1 


0.54  +  0.46  cos(irf/T);  -t/2  <  f  <  t/2 
0  ;  elsewhere 


(4.19) 


Note  that  the  Hamming  window  is  symmetric,  therefore,  when  convolved 
with  a  waveform,  the  result  is  equivalent  to  shifting  the  window  to 
specific  locations  of  the  waveform,  then  taking  the  weighted  average  of 
neighboring  points  covered  under  the  window,  as  shown  in  Fig.  4.5.  The 
complex  data  points  can  be  interpolated  and  smoothed  by  convolving  the 
Hamming  Window  with  the  data  itself  as  shown  in  the  following  equation. 


c(fj> 


m 

£  H(fkjfjtT)  Ws( fk)c( fk) 
k=n _ 

"  H(fk,fj,T)  Ws(fk) 
k=n 


(4.20) 


where 


fj  —  t/2  <  fn s  .  .  .  ,  fjjj  <  fj  +  t/ 2j  m  >  n 
H(fk,fj,x)  =  0.46  +  0.54  cos(1T[fk~fj)) 

T 

Ws(fk)  =  (fk+l  "  fk-l)/2  (4.21) 
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The  value  of  the  waveform  at  f  =  fj,  i.e.,  c(fj)  is  estimated  by 
calculating  the  weighted  average  of  the  data  values  covered  under  the 
Hamming  window  Hff^.fj.t).  The  weight  for  each  data  point  under  the 
window  is  determined  by  the  Hamming  window;  the  further  the  data  from 
the  center  the  smaller  the  weight. 

Note  that  the  weighting  function  Ws(f|<)  is  added  to  compensate  for 
possible  errors  caused  by  uneven  distribution  of  data  points,  as  will  be 
shown  numerically  in  an  example.  The  term  in  the  numerator  of  Eq. 

(4.20)  is  a  normalization  factor  so  that  the  sum  of  the  weights  is  1. 

The  necessity  of  using  the  normalization  factor  becomes  evident  when  the 
waveform  c(fx)  is  constant.  In  this  case  the  interpolated  values  c(fj) 
remain  constant,  as  expected. 

Now  let's  interpret  Eq.  (4.20)  in  the  time  domain.  The  inverse 
Fourier  Transform  of  the  Hamming  Window  is  a  low-pass  filter  given  by 
[35] 


H(t,x)  =  Tsinftm:)  .  (0.54  -  0.08(ti)2) 

_T_ - - —? - 

1  -  (  t  X ) 


(4.22) 


The  convolution  in  Eq.  (4.20)  therefore  is  equivalent  to  a  time-domain 
low  pass  filtering  as  shown  in  Fig.  4.6.  Note  that  the  time-domain 
Hamming  window  is  real  and  its  amplitude  has  a  maximum  at  t=0,  and  is 
about  8  db  down  at  t  -  1/t.  Since  our  calibration  scheme  as  shwon  in 
Eq.  (4.2)  sets  the  phase  center  to  be  located  at  the  center  of  the 
targets,  the  time-domain  returns  at  time  equal  to  zero  would 


Fig.  4.6.  Interpretation  of  the  Interpolation  in  the  time  domain. 

The  Hemming  window  becomes  a  low-pass  filter  in  the  time 


correspond  to  the  centers  of  the  targets.  The  time  t  and  the  target 
range  R  are  related  by 

t  =  2R/c  (4.23) 

where  c  is  the  speed  of  light.  For  a  target  of  length  i ,  a  major  part 
of  the  time-domain  returns  is  confined  to  a  zone  given  by 

-  i/c  <  t  <  i/c  (4.24) 

In  addition  most  of  the  clutter  is  greatly  removed  in  distance  (or  time) 
from  the  target.  Mote  that  the  clutter  here  is  very  different  from 
random  noise,  and  is  confined  to  specific  mechanisms  such  as  reflections 
from  walls.  Therefore,  the  clutter  can  be  removed  by  the  Hamming  window 
if  the  window  width  is  chosen  properly.  For  example,  if  x  is  chosen  as 
300  MHZ,  then  1/r  is  equivalent  to  a  range  of  50  cm.  Since  model 
lengths  are  about  10  cm,  a  major  part  of  the  returns  comes  from  the  zone 
below  (see  Fig.  4.6). 

-  0.05  I  <  t  <  0.05  I  (4.25) 

T  T 

Therefore  most  clutter  will  br  removed.  The  effectiveness  of  this 
algorithm  is  illustrated  in  the  following  example. 

Suppose  we  are  given  6  data  points,  which  are  sampled  from  a 
waveform  defined  as 

c(f)  =  f 
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The  waveform  is  a  ramp,  as  shown  in  Fig.  4.7.  Note  that  the  dots  on  the 
the  waveform  indicate  sample  points.  Note  that  the  spacing  between  each 
sample  is  not  constant.  Suppose  that  we  want  to  estimate  the  value  of 
the  waveform  at  f  =  3  based  on  the  6  data  points.  The  true  value  at 
f  =  3  should  be  C  =  3.  Now  we  apply  a  Hamming  window  on  the  data  points 
with  its  center  located  at  f  =  3.  Since  there  are  four  data  points  on 
one  side  of  the  window,  and  only  two  on  the  other  side,  the  window  is 
unbalanced,  i.e.,  the  contribution  from  one  side  is  larger  than  the 
other.  From  Eq.  (4.20),  the  waveform  value  at  C  =  3  can  be  estimated  as 


c(3)  =  2.9489 


and  the  error  is  given  as 


%  error  =  c(3)  -  c(3)  =  1.7% 


Now  suppose  the  weighting  function  Wjffy)  is  not  used,  despite  the  data 
points  are  not  equally  spaced,  then 


c ( 3 )  =  3.35 


and  the  error  is  given  as 


%  error  =  c ( 3 )  -  c  ( 3 )  =  11.7% 
cTTJ 


which  is  several  times  higher.  The  weighting  function  Ws(f|<)  in 
Eq.  (4.20)  is  obviously  very  effective. 


■  ■  .  '  -  1  .  V  -  yi n'.V.v.s'f.An  * 


Fig.  4.7.  An  example  of  the  interpolation.  There  are  6  data  points 
located  at  f  =  1,  3,  3.5,  4,  4.5,  and  5.  The  underlying 
function  is  c(f)  =  f.  It  is  desired  to  estimate  the  data 


value  at  f  =  3. 


E.  POST-PROCESSED  RADAR  RETURNS 


The  accuracy  of  the  post-processed  data  can  be  examined  by 
comparing  processed  returns  of  spheres  with  the  theoretical  returns  of 
the  spheres.  An  example  is  given  in  Fig.  4.8  for  a  3.5  inch  (diameter) 
sphere  calibrated  by  a  6  inch  reference  sphere.  It  is  seen  that  the 
error  in  the  amplitude  waveforms  is  small.  The  error  in  the  phase  wave¬ 
forms  is  also  small  except  for  a  ramp-type  phase  offset.  The  offset 
is  caused  by  a  difference  in  positions  between  the  center  of  the  3.5 
inch  sphere  and  the  6  inch  sphere.  Based  on  this  type  of  comparison,  we 
found  that,  in  general,  the  accuracy  of  the  measured  data  is  enough  for 
our  application  here. 

Plots  of  the  processed  radar  returns  of  five  model  ships  and  four 
model  aircraft  are  shown  in  Appendix  A.  The  aspect  angles  of  the 
targets  are  0,  90  and  180  degrees,  and  the  polarization  is  vertical. 
Smoothed  and  unsmoothed  waveforms  are  shown  at  the  same  time  for  compar¬ 
ison.  Note  that,  for  this  modelling,  the  amplitude,  frequency,  and 
phase  are  scaled  as  follows: 

As  =  An 
fs  =  f/n 
es  =  9 

where 

n  =  scale  factor 
As  =  scaled  amplitude 
fs  =  scaled  frequency 
0S  =  scaled  phase. 
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Fig.  4.8.  Comparison  between  measured  and  theoretical  returns  of  a 


3.5  inch  sphere  at  f  =  2-4GHZ 


\% 


E 


CHAPTER  V 


TESTING  OF  CLASSIFICATION  ALGORITHMS 


A.  INTRODUCTION 


This  chapter  will  discuss  the  experiment  which  were  conducted  to 
test  the  performance  of  the  target  classification  algorithms  developed 
in  Chapter  III. 

Figure  5.1  is  a  flowchart  showing  the  experimental  process  for 
testing  the  classification  algorithms.  As  shown  in  the  figure  a  data 
base  is  comprised  of  mul ti -frequency  amplitude  data  A^ ( f j )  and  phase 
data  QiJfj).  where  K  means  the  Kth  target,  and  fj  means  the  jth 
frequency.  Two  data  bases  were  generated,  each  consisting  of  multi  - 
frequency  radar  data  for  five  ships  and  four  aircraft.  Techniques  for 
generating  these  two  experimental  data  bases  have  been  discussed  in 
Chapter  VI,  section  C.  Based  on  the  data  bases,  we  generated  catalogs 
of  features  for  time-domain  classification  and  for  frequency-domain 
classification.  Features  for  test  targets  were  generated  by  injecting 
noise  into  the  data  stored  in  the  data  base.  In  section  B  we  will 
discuss  the  generation  of  the  test  data  and  the  noise  model  used.  The 
test  target  data  were  then  classified  either  in  the  frequency  domain  or 
in  the  time  domain,  as  discussed  in  Chapter  III.  By  performing  the 
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Figure  5.1. 


A  flowchart  showing  the  experimental  process  for 
testing  the  classification  algorithms. 


classifications  of  test  target  data  several  times,  where  each  experiment 
Is  statistically  Independent,  we  can  obtain  an  estimate  of  the 
probability  of  mlsclasslflcatlon.  Section  C  will  discuss  the 
statistical  techniques  for  estimating  the  probability  of 
mlsclasslflcation.  In  section  D  we  will  present  experimental  results 
for  the  classification  of  ships,  and.  In  section  E,  results  for  the 
classification  of  aircraft. 

B.  NOISE  MODEL 

There  are  a  lot  of  factors  which  affect  the  noise  present  In  the 
post-processed  radar  returns.  In  general,  noise  models  which  can 
precisely  describe  the  characteristics  of  the  noise  In  practice  are 
too  complex  to  be  treated  with  available  mathematical  methods,  or  are 
too  cumbersome  to  compute. 

Consider  three  noise  models  as  shown  below. 

1.  s  =  s  +  e0  (additive  noise) 

2.  s  =  ejs  (multiplicative  noise) 

o*  /V 

3.  s  =  e2S  +  £3  (additive  and  multiplicative  noises) 

where,  s  is  the  original  signal,  s  Is  the  noise-contaminated  signal,  and 

/v  >v  a# 

eg,  ej,  £2  and  £3  are  Gaussian-dlstrlbuted  random  noise  terms. Based  on 
the  noise  models  and  a  data  base,  we  calculated  the  probabilities  of 
mlsclasslflcatlon  as  a  function  of  signal -to-noise  power  ratio.  The 
result  showed  that  the  probability  of  mlsclasslflcatlon  for  each  model 
for  a  given  signal -to-noise  ratio  was  similar  to  each  other.  This 
Indicates  that  the  total  power  of  the  noise,  rather  than  Its  specific 


r 

v 


characteristics,  has  more  Impact  on  the  performance  of  the 
classification  algorithms.  Furthermore,  since  there  are  numerous 
Independent  noise  sources,  none  of  which  is  dominant,  the  statistics  of 
the  sum  of  these  noises  can  be  approximately  described  as  Gaussian  by 
virtue  of  the  Central  Limit  Theorem.  Therefore  an  additive  Gaussian 
noise  model  would  be  appropriate  for  our  application  here. 

A  noise-free  data  point  at  a  particular  frequency  is  a  vector 
consisting  of  an  in-phase  component  I  and  a  quadrature-phase  component 
Q,  which  are  defined  as 

I  =  A  cos  0 

Q  =  A  sin  0  (5.1) 

where  A  is  the  amplitude  and  0  is  the  phase.  The  in-phase 
component  and  the  quadrature-phase  component  of  a  noise-contaminated 
data  point  can  be  generated  by  adding  Gaussian  noise  to  the  noise-free 
data  point,  i.e., 

?  =  I  +  nj 

Q  =  Q  +  nQ  (5.2) 

where  nj  and  nQ  are  independent,  Gaussian-dlstributed  random  variables 
with  zero  mean  and  variance  c^.  The  relation  between  the  noisy  data 
point  and  the  noise-free  data  point  is  Illustrated  in  Figure  5.2,  where, 
on  a  I-Q  plane,  It  is  shown  that  the  noise-free  data  point  may  be 
scattered  by  the  noise  vector  over  a  circular  region  (shaded).  Here  we 
assume  that  the  power  of  each  noise  components,  nj  and  ng,  is  the  same 
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so  that  the  shaded  region  Is  circular.  It  Is  noted  that  the 
amplitude  of  the  noise  vector  Is  given  as 


An  =  J nx2  +  n* 


(5.3) 


which  Is  Rayleigh-distributed,  and  the  phase  of  the  noise  vector  Is 
given  as 


en  =  tan'1  ® 


(5.4) 


which  Is  uniformly  distributed  [31].  From  Equation  (5.2)  we  can 
calculate  a  noisy  amplitude  and  a  noisy  phase  as 


/~2  ~2 
A  =  -\/  I  +  Q 


(I  +  nj)2  +  (Q  +  nQ)2 


(5.5) 


9  =  tan'1  (i] 


=  tan'1^  } 


(5.6) 


The  noisy  phase  in  the  above  equation  can  be  used  to  calculated  a  noisy 
W,  l.e.. 


•G? 


-1  mJi  aJ  aJCjkJi  fc-.4!  JL.*  •.  »_**.  «.  1* fc_*\  A  **.V  *. 


Wi  =  ei*i  -  ei+lxl+l 


(5.7) 


In  the  case  where  only  amplitude  Is  available  the  amplitude  in  the 
presence  of  the  noise  is  given  by 

(5.8) 

Noisy  time-domain  data  can  be  obtained  by  taking  the  inverse  discrete 
Fourier  transform  of  the  noisy  frequency-domain  data  defined  in  Eq. 
(5.2).  ^ 

All  of  the  noise-contaminated  features  above  are  defined  based  on 
the  additive  Gaussian  noise  model  defined  in  Eq.  (5.2).  A 
signal-to-noise  power  ratio  (SNR)  for  such  a  model  can  be  defined  as  the 
ratio  of  total  signal  power  and  total  noise  variance,  i.e., 

S  =  I2  +  Q2 
¥  VAR(ni)  +  VAR(nq) 

=  A2  (5.9) 

2 

2a 

In  our  classification  experiments,  the  signal  power  A2  was  first 
estimated  using  the  data  in  the  data  base.  By  varying  the  noise  power 
a2  the  SNR  was  varied  over  a  range,  and  the  probability  of 
mi sclassifi cation  was  determined  for  each  SNR  value. 


C.  ESTIMATION  OF  THE  PROBABILITY  OF  MISCLASSIFICATION 

Analytical  calculations  of  error  probabilities  for  our 
classification  schemes  here  are  rather  involved  because  they  may  require 
integrations  of  multi -dimensional  probability  density  functions  with 
extremely  complicated  boundaries.  Furthermore  the  density  functions  may 
not  be  readily  obtainable.  Therefore  one  has  to  resort  to  statistical 
methods  for  estimating  the  probability  of  mi sclassifi cation. 

Suppose  there  are  M  targets,  and  the  probability  of  occurrence  of 
each  target  Cj,  i.e.,  the  a  priori  probability,  is  given  by  p(cj), 
j  =  1,  2,  .  .  .  ,  M.  The  probability  of  mi sclassifi cation  can  be  given 
as 

M 

P(Error)  =  £  P(Error/cj)P(c-s)  (5.17) 

J-l 

where  P(Error/cj)  is  the  conditional  probability  of  error  given  target 
j.  A  good  statistical  estimate  for  the  conditional  probability 
p(Error/cj)  is  the  Maximum  Likelihood  Estimate  (MLE)  [24]. 

Now  let's  derive  the  MLE  of  P(Error/cj).  Let 

P(Error/cj)  h  Pj  (5.18) 

Suppose  nj  independent  classification  experiments  are  performed  on 
target  j.  This  can  be  done  by  using  different  seeds  for  generating 
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noises  In  each  experiment  so  that  the  noises  are  Independent  of  each 
other.  In  each  experiment  the  probability  of  mlsclassifying  target  j  is 
Pj.  It  is  evident  that  the  probability  of  misclasslfying  target  j  mj 
times  can  be  described  by  a  binomial  probability  distribution,  i.e.,  mj 
is  Binomial -distributed  with  parameter  pj,  and 


/  /  i  "Vj  nj"mj 

p{mj/pj}  =  (m  )  p  (1-pj) 

J 


(5.19) 


A  feasible  choice  of  the  Likelihood  Function  would  be  the  term  p{ mj /pj } 
defined  by  Eq.  (5.19)  [24].  The  MLE  of  pj  is  the  one  which  maximizes 
the  Likelihood  Function,  i.e.. 


or 


Pj  =  max  ptmj/pj} 


3  log  p  (mj/Pj) 


(5.20) 


Pj=Pj 


-  0 


(5.21) 


where  pj  is  the  MLE  of  pj.  Now 

3  [log  p  (mj/pj}] 

3Pj 


"j, 


m, 


=  3  log  {(  “)  p.  (1-pj)  “} 

IpJ  mj  j 

nj 

=  3  (log(  )  +  mjlogpj  +  (nj-mj)log(l-pj)} 

=  !!i  - 
Pj  (1-PjT 
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(5.22) 


J.  .■j.'-V..1.  £  -  .. 
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From  Eq.  (5.21)  and  Eq.  (5.22)  we  have 


i 


m, 


n 


j 


nj~mj  =  0 
d-Pj) 


or 

pj  =  "J  (5 

"j 

When  the  sample  size  nj  is  large  (normally  njpj  >  5),  a  (1  -  a) 
confidence  interval  for  pj  can  be  approximately  given  as  [24] 

Pj  ±  5a-\/v AR(pj) 


or 


Pj  <  Pj 


a) 


where  Za  is  defined  by  $Ua)  =  l-ot/2 ,  where  t  is  the  standard  normal 
distribution,  and 


*(5a)  = 
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If 


?E  1  Or^h) 

where  g  a  is  defined  as  *Ua)  =  l-a/2. 


Now 

VAR(Pe) 


M  - 

VAR  (1  I  p-t) 
ffj=l 

M 

1  Z  VAR(p0) 
"2  j=l 
M 


=  1  VAR{p0) 

w 


=  1  Pod-Po^ 
FT  n 


A  A 


where  p^,  P2,  ....  pm  are  assumed  to  be  mutually  independent  and 


VAR ( p i )  =  VAR(p2)  =  .  .  .  =  VAR(pm)  =  VAR ( p0 ) 


Therefore  the  (1-a)  confidence  interval  for  Pe  becomes 
Pf 


£* 


’E  ±  g^/PQ^-Po) 


(5.25) 


An  algorithm  for  estimating  the  probability  of  misclassification  is 
given  below. 

Suppose  there  are  M  data  files  in  the  data  base  representing, 
respectively,  M  classes,  and  each  containing  N  complex  frequency 
components,  i.e.,  the  I-data  and  Q-data. 
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1.  Create  M  targets  by  injecting  noise  into  the  M  data 
files  respectively,  i.e.,  for  each  I  data  and  Q  data 
of  each  of  the  N  frequencies  of  each  file,  generate 
two  Gaussian  random  numbers,  and  contaminate  the  data 
as  shown  in  Eq.  (5.2).  Use  different  seeds  for  each 
file. 

2.  Classify  the  M  targets  according  to  the  algorithms  in 
Section  C  of  Chapter  III  (frequency-domain  classification) 
and  in  Section  D  of  Chapter  III  (time-domain 
classification). 

3.  Repeat  the  above  experiments  (steps  1  and  2)  a  number 
of  times,  using  different  seeds  for  each  experiment. 

4.  Calculate  the  probability  of  misclassification  using 
Eq.  (5.24). 

D.  CLASSIFICATION  OF  SHIPS 

In  this  section  we  will  discuss  our  experiments  on  the 
classification  of  ships  and  the  results. 

The  experiments  made  use  of  a  data  base  consisting  of  15  classes  as 
defined  in  Table  5.1.  The  aspect  angles  of  the  ships  were  chosen 
to  be  0,  90,  and  180  degrees  so  that  data  of  the  same  ships  at 
different  aspect  angles  (e.g.,  class  1,  2,  and  3)  would  be  independent 
of  each  other.  Since  data  of  different  ships  (e.q.,  class  1  and  4) 
were  also  independent,  all  classes  in  the  data  base  would  be 
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independent  of  each  other.  The  number  of  frequencies  N 
employed  for  the  classifications  were  N  =  2,  4,  and  8,  with  a 
frequency  increment  of  Af  =  0.4  MHZ.  Young  and  Walton  found  [13]  that 
the  duration  of  a  time-domain  ramp  response  R(t)  was  less  than  6 
transit  time  across  the  length  i  of  a  target,  i.e., 

R(t,t0)  »  0  for  t  <  t0  or  t  >  t0  +  (6fc/c)  (5.26) 

where  C  is  the  speed  of  light  and  t0  is  a  time  reference.  For  the  above 
time-limited  function,  the  Shannon's  sampling  theorem  requires  that  the 
frequency-sampling  rate  Af  satisfy 

Af  <  c/6 l  (5.27) 

The  ships  have  dimensions  of  the  order  of  100  meters,  which  implies  that 
the  sampling  rate  for  the  ship  returns  should  satisfy 

Af  <  0.5  MHZ  (5.28) 

Based  on  the  criterion  in  Eq.  (5.28)  we  chose  the  frequency-increment  Af 
as  0.4  MHZ.  The  minimum  frequency  available  in  the  data  base  was  1.67 
MHZ,  therefore,  the  frequency  ranges  were  chosen  as  1.7  MHZ  <  f  <  2.1  MHZ 
for  N  =  2,  1.7  MHZ  <  f  <  2.9  MHZ  for  N  =  4,  and  1.7  MHZ  <  f  <  4.5  MHZ 
for  N  =  8.  The  corresponding  range  of  wavelengths  for 
1.7  MHZ  <  f  <  4.5  MHZ  is  67  m  <  A  <  176  m,  and  the  ratios  of  the  ship's 
dimensions  D  and  the  wavelengths  are  of  the  order  of 


D/A  -  1 


(5.29) 


which  means  the  radar  returns  are  in  the  lower  Resonance  region.  For 
frequency-domain  classifications  the  parameter  3  in  Eq.  (3.15)  was 
chosen  as  3  =  »,  3=0,  and  3=1,  corresponding,  respectively,  to  uses 
of  'amplitude  features  only',  ‘phase  features  only',  and  'both  amplitude 
and  phase  features'.  A  feature  called  'relative  amplitude'  was  also 
employed  for  classifications  in  the  frequency  domain,  and  was  defined  as 
follows. 

Xj  =  A-j/Ai+i  ;  i  =  1,  2,  .  .  .  ,  N  -  1  (5.30) 

where 

Xj  =  ith  relative  amplitude 

Aj  =  ith  amplitude 

N  =  number  of  frequencies 

The  division  undertaken  in  Eq.  (5.30)  was  designed  for  the  removal  of 
possible  multiplicative  bias  embedded  in  the  amplitude  data. 

For  time-domain  classification,  the  Discrete  Inverse  Fourier 
Transform  (DIFT)  was  applied  to  obtain  real  numbers  for  the  time-domain 
signals.  Frequency-domain  data  were  assumed  to  be  symmetric  around 
f  =  0,  i.e. , 
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c(f)  =  c*(-f) 


(5.31) 


where,  C*(f)  Is  the  complex  conjugate  of  the  frequency-domain  data  C(f), 
so  that  the  corresponding  time-domain  signals  would  be  real  numbers. 

It  Is  known  that  we  can  reduce  the  error  probability  If  the  target 
aspect  angles  are  known.  In  practice  the  aspect  angles  should  be 
estimated  with  some  accuracy  otherwise  the  probability  of 
ml sclasslfl cation  may  Increase  sharply.  Let's  consider  the  effects  of 
an  error  of  15  degrees  in  the  estimated  target  aspect  angles,  as  shown 
In  Fig.  5.3.  The  change  In  path  length  aP  for  a  change  In  aspect  angles 
by  15  degrees  is  given  by 

Ap  =  2  •  AJl 

=  2  •  (2.  -  £  cos  15°) 

2  2 

or 


Ap  =  3.4%  D  (5.32) 

which  Is  only  a  fraction  of  the  total  target  dimension  D.  Consider  a 
data  base  as  shown  In  Table  5.2,  where  the  aspect  angles  are  15  degrees 
apart.  Tables  5.3  -  5.5  are  matrices  which  show  the  N.N.  distances 
between  any  two  classes  In  the  data  base  for  6  =  »  (amplitude  only). 
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TABLE  5.3 


A  MATRIX  SHOWING  THE  DISTANCE  BETWEEN  ANY  TWO  CLASSES. 

THE  ENTRY  AT  ROW  I  AND  COLUMN  J  IS  THE  N.N.  DISTANCE  BETWEEN  CLASS  I 
AND  CLASS  J.  ONLY  AMPLITUDE  DATA  ARE  USED  TO  CALCULATE  THE  DISTANCE. 


CLAS 

S  1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

II 

12 

13 

14 

15 

1 

0. 

18. 

23. 

162. 

155. 

150. 

276. 

291. 

203. 

62. 

53. 

20. 

166. 

181. 

230. 

2 

18. 

0. 

10. 

164. 

158. 

157. 

277. 

294. 

205. 

69. 

57. 

30. 

164. 

175. 

232. 

3 

23. 

10. 

0. 

166. 

160. 

160. 

280. 

297. 

206. 

66. 

54. 

32. 

167. 

175. 

231. 

4 

162. 

164. 

166. 

0. 

20. 

89. 

199. 

206. 

125. 

193. 

191. 

179. 

120. 

140. 

15B. 

5 

155. 

158. 

160. 

20. 

0. 

73. 

199. 

204. 

121. 

186. 

1 84. 

170. 

119. 

144. 

153. 

6 

150. 

157. 

160. 

89. 

73. 

H 

225. 

205. 

103. 

191. 

1B9. 

160. 

109. 

152. 

147. 

7 

276. 

277. 

280. 

199. 

199. 

225. 

0. 

105. 

173. 

288. 

286. 

284. 

183. 

192. 

218. 

8 

291. 

294. 

297. 

206. 

204. 

205. 

105. 

0. 

147. 

312. 

31 D. 

297. 

176. 

199. 

204. 

9 

205. 

208. 

125. 

121. 

103. 

173. 

147. 

0. 

23  8. 

235. 

212. 

65. 

107. 

146. 

10 

62. 

69. 

66. 

193. 

186. 

191. 

288. 

312. 

23  8. 

0. 

15. 

57. 

206. 

213. 

25B. 

II 

53. 

57. 

54. 

191. 

184. 

189. 

286. 

310. 

235. 

15. 

0. 

48. 

200. 

207. 

257. 

12 

20. 

30. 

32. 

179. 

170. 

160. 

284. 

297. 

212. 

57. 

48. 

0. 

175. 

190. 

236. 

13 

166. 

164. 

167. 

120. 

119. 

109. 

183. 

176. 

65. 

206. 

200. 

175. 

0. 

55. 

134. 

14 

181. 

175. 

175. 

140. 

144. 

152. 

192. 

199. 

107. 

213. 

207. 

190. 

55. 

0. 

152. 

15 

230. 

232. 

231. 

158. 

153. 

147. 

218. 

204. 

146. 

258. 

257. 

236. 

134. 

152. 

0. 

218.  204.  146 


258.  257.  236 


134.  152 


0 


8=0  (phase  only),  and  8=1  (amplitude  and  phase).  Table  5.6  is  a 
matrix  which  shows  the  cross-correlation  coefficients  of  any  two  classes 
in  the  data  base.  Note  that  the  N.N.  distance  and  the  cross-correlation 
coefficient  between  a  class  and  itself  are  0  and  1,  respectively,  as 
shown  in  the  diagonal  of  the  matrices.  In  the  tables  the  nearest 
neighbor  of  a  ship  other  than  the  ship  itself  is  underlined  once  if  the 
nearest  neighbor  comes  from  the  same  ship,  and  is  underlined  twice  if  it 
comes  from  different  ships.  The  percentage  of  entries  with  double 
underlines  in  the  matrices  indicates  the  classification  error  introduced 
by  an  error  of  15  degrees  in  the  estimated  target  aspect  angles.  The 
percentage  of  double  underlines  is  13%  in  Table  5.3,  and  20%  in  Tables 
5.4,  5.5,  and  5.6.  This  indicates  the  tolerable  uncertainty  of  target 
aspect  angles.  If  the  mi sclassifi cation  probability  is  too  high  the 
target  aspect  angles  should  be  determined  more  accurately. 

The  signal  power  of  the  data  in  the  data  base  may  vary  depending  on 
frequency,  aspect  angles,  and  geometries  of  the  targets.  The  noise 
added  to  each  data  point  was  set  to  have  the  same  variance  so  that  data 
with  higher  signal  powers  would  have  higher  signal-to-noise  power 
ratios.  The  signal  power  of  each  class  was  estimated  by  taking  the 
average  over  all  data  points  of  that  class,  i.e., 

N 

p  =  1  Z  ( f ^ ) ;  A  =  Amplitude  (5.33) 

¥  i=l 

Since  additive  zero-mean  Gaussian  noise  is  used,  the  noise  power  is 
equal  to  its  variance  a?.  The  average  signal-to-noise  ratio 


TABLE  5.6 

A  MATRIX  SHOWING  THE  CROSS-CORRELATION  COEFFICIENT 
FOR  ANY  TWO  CLASSES 


CLASS  I  2 


3  4 


6  7 


9  10 


13  14  15 


1 

1.00  0.  81  0.75 

0.58  0.59  0.44 

0.61  0.51  0.49 

0.68  0.75  0.88 

0.51  0.71  0.57 

2 

0.81  1.00  0.90 

0.66  0.5B  0.52 

0.76  0.54  0.83 

0.70  0.77  0.87 

0.79  0.54  0.55 

3 

0.75  0.90  1 .00 

0.83  0.7  8  0.55 

0.66  0.47  0.51 

0.B6  0.80  0.78 

0.63  0.59  0.56 

4 

0.58  0.66  0.83 

1 .00  0.9B  0.55 

0.54  0.42  0.50 

0.69  0.66  0.5B 

0.55  0.55  0.55 

5 

0.59  0.58  0.78 

0.98  1.00  0.61 

0.47  0.48  0.57 

0.67  0.65  D.57 

0.56  0.57  0.57 

6 

0.44  0.52  0.55 

0.55  0.61  1.00 

0.61  0.72  0_173 

0.55  0.53  0.48 

0.67  0.50  0.49 

7 

0.61  0.76  0.66 

0.54  0.47  0.61 

1.00  0.70  0.77 

0.50  0.57  0.62 

0.88  0.68  0.66 

8 

0.51  0.54  0.47 

0.42  0.48  0.72 

0.70  1.00  0^91 

0.45  0.41  0.61 

0.70  0.50  0.63 

9 

0.49  0.63  0.51 

0.50  0.57  0.73 

0.77  OJM  1.00 

D.50  0.52  0.44 

0.82  0.65  0.70 

10 

0.68  0.70  0.86 

0.69  0.67  0.55 

0.50  0.45  0.50 

1.00  0.98  0.79 

0.61  0.61  0.53 

II 

0.75  0.77  0.80 

0.66  0.65  0.53 

0.57  0.41  0.52 

0,98  1.00  0.89 

0.65  0.62  0.48 

12 

0.88  0.87  O.’B 

0.5B  0.57  0.48 

0.62  0.61  0.44 

0.79  0.89  1.00 

D.66  0.67  0.42 

13 

0.51  0.79  0.63 

0.55  0.56  0.67 

0.88  0.70  0.82 

0.61  0.65  D.66 

1.00  0.62  D.72 

14 

0.71  0.54  0.59 

0.55  0.57  0.50 

0.68  0.50  0.65 

0.81  0.62  0.67 

0.62  1.00  0.B1 

15 

0.57  0.55  0.56 

0.55  0.57  0.49 

0.66  0.63  0.70 

0.53  0.48  0.42 

0.72  0.81  1.00 

»»:  -.•••  v.-v  .  ■ 


is  given  by 


N 

SNR  =  1  E  (Pi/a2) 

IT  j=l 

k 

=  1_  ( 1  Z  Pj)  (5.34 

2  S’  J=1 

a 

where  pj  is  the  average  power  of  class  j,  and  K  is  the  number  of 
classes,  which  is  equal  to  15  here.  The  term  in  the  parenthesis  in  Eq. 
(5.34)  is  the  average  power  over  all  classes  in  the  data  base,  and  was 
calculated  to  be  approximately  44  dB-square-meters.  The  probability  of 
mi sclassifi cation  was  calculated  based  on  Equation  (5.24)  with  the 
number  of  experiments  =  15. 

Let's  use  Equation  (5.25)  to  calculate  a  90  percent  confidence 
interval  for  the  probability  of  mi  sclassifi  cation  Pfr.  Let  the 
probability  of  misclassification  be  about  30  percent,  i.e.,  p0  ~.3 
With  M=15,  n=15,  a=0.1,  and  ?a  satisfying 

<&(?)  =  1  -  «  =  0.95  . 

«  7 

or 

=  1.65  . 


the  confidence  interval  for  Pg  is  given  as 


P 


E 


+ 


5  /Pod-Po) 
a  V - Wn - 


=  P 


+ 

E  " 


1.65 


VO. 3(1-0. 3) 
17x175 - 


LO 


=  P£  ±  5  percent  . 


(5.35) 


Therefore  a  90  percent  confidence  interval  for  Pf  «  30  percent  is  ±  5 
percent. 

In  the  following  we  will  present  plots  of  mi sclassifi cation 
percentages  as  a  function  of  signal-to-noise  power  ratio  (S/N)  under  the 
following  conditions: 

1.  As  a  function  of  the  number  of  frequencies  utilized, 

2.  For  different  classification  algorithms, 

3.  With  or  without  a  priori  knowledge  of  target  aspect  angles. 

Figures  5. 4-5. 7  are  plots  of  error  probabilities  with  the  number  of 
frequencies  N  as  a  parameter.  In  all  four  figures  we  assume  that  target 
aspect  angles  are  known  a  priori,  and  the  number  of  frequencies  N  is  set 
to  2,  4  and  8.  Figures  5.4,  5.5,  and  5.6  are  results  of  frequency- 
domain  classifications  using,  respectively,  amplitude  returns  only, 
phase  returns  only,  and  both  amplitude  and  phase  returns.  Figure  5.7 
shows  the  result  for  time-domain  algorithms.  It  can  be  seen  that  the 
error  probabilities  can  be  improved  significantly  by  increasing  the 
number  of  frequencies.  This  is  especially  true  for  time-domain 
classi fications. 

Figures  5.8-5.10  show  the  effects  of  employing  different  types  of 
classification  algorithms.  For  convenience  let  A  stand  for  "use  only 
amplitude",  W  for  "use  only  phase",  A  and  W  for  "use  both  amplitude  and 
phase",  and  TIME  for  "use  time-domain  algorithms".  In  all  three  figures 
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-12.  -9.  -6.  -3.  0.  3.  6.  9.  12.  15.  18.  21.  29.  27. 

POST-PROCESSING  S/N  IN  D3 

Figure  5.4.  Plots  of  error  probabilities  with  the  number  of 

frequencies  N  as  a  parameters.  Only  amplitude  is  used 
Aspect  angles  are  assumed  known.  A  90%  confidence 
interval  for  a  classification  error  of  30%  is 


apDroximately  ±  5%  . 


we  assume  that  target  aspect  angles  are  known  a  priori.  The  number  of 
frequencies  used  is  N  =  8  in  Figure  5.8,  N  =  4  in  Figure  5.9,  and  N  =  2 
in  Figure  5.10.  It  can  be  seen  that  the  time-domain  algorithms  produce 
the  best  results  in  Figure  5.8  and  Figure  5.9.  In  Figure  5.10  where 
N  =  2,  the  time-domain  results  are  worse  than  those  of  'A'  or  'A  or  W  . 

It  can  be  seen  from  the  definition  of  the  noise  model  that  certain 
amounts  of  noise  introduced  in  the  I  and  Q  components  may  cause  large 
errors  in  VI  due,  in  part,  to  branch-cut  effects.  As  a  result,  the  error 
in  the  phase  features  dominates  the  decision.  The  error  becomes  larger 
if  the  number  of  frequencies  is  increased.  Due  to  this  effect,  it 
can  be  seen  in  Figure  5.8,  where  N  =  8,  the  'A'  curve  is 
better  than  'A  and  W  or  'W  curves.  In  Figure  5.10  where 
the  number  of  frequencies  becomes  smaller  (N  =  2)  and  the  'A  and  W 
curve  becomes  better  than  the  'A1  curve. 

Figure  5.11  shows  the  error  probabilities,  in  which  relative 
amplitudes,  "A,  as  defined  in  Equation  (5.30),  are  employed.  As  is  the 
case  for  the  W’s,  the  relative  amplitudes  may  also  be  distorted 
significantly  when  certain  amounts  of  noise  are  introduced  into  the  I 
and  Q  components,  resulting  in  the  same  effects  as  discussed  above.  In 
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Figure  5.11  we  can  see  that  the  ' W 1  curve  is  better  than  the  'A  and  W 
curve  which  indicates  that  errors  in  the  "ft's  are  even  larger  than  those 
in  the  W's.  A  conclusion  we  can  draw  here  is  that,  when  using 
incompatible  features  such  as  the  Vs  and  W's,  we  must  carefully  select 
the  parameters  so  that  the  features  have  equal  effects  in  making  the 
N.N.  decision. 

Figures  5.12-5.15  are  comparisons  of  error  probabilities  between 
cases  where  target  aspect  angles  are  either  known  or  unknown. 

Figures  5.12-5.15  show,  respectively,  the  results  of  using  amplitude 
returns  only,  phase  returns  only,  both  amplitude  and  phase  returns,  and 
the  time-domain  algorithms.  Note  that  if  the  aspect  angles  are  known 
the  maximum  error  is  4/5  =  80  percent,  and  if  the  aspect  angle  are 
unknown  the  maximum  error  is  14/15  =  93  percent.  It  is  seen  that 
knowledge  of  target  aspect  angles  will  decrease  the  errors. 

Based  on  the  above  figures  we  can  derive  a  list  of  minimum  S/N 
values  required  for  10  percent  or  less  error  probabilities  as  shown  in 
Table  5.7.  This  gives  us  a  quantitative  indication  of  the  performance 
of  each  algorithm.  From  the  table  we  note  that,  for  N  =  2,  the  'A  ana  W 
algorithm  requires  the  least  S/N,  followed  by  'A  '  ,  'TIME1,  and 

' W 1 ,  in  increasing  order.  For  N  =  4,  the  required  SNR  are  about  the 
same  for  the  'A  and  W' ,  'A',  and  'TIME*  algorithms,  while  * W '  algorithms 
requires  higher  S/N.  For  N  =  8,  the  minimum  SNR  required  for  each 
algorithm  has  the  following  order: 

'TIME'  <  'A'  <  'A  and  W'  <  '£  &  W  »  ' W  < 
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TABLE  5.7 

MINIMUM  S/N  REQUIRED  (in  dB)  FOR  10  PERCENT  OR  LESS 
CLASSIFICATION  ERRORS.  FOR  EACH  ALGORITHM  (e.g.,  A&W) 
THERE  ARE  FOUR  S/N  VALUES  (e.g.,  19,  17,  13,  15) 
CORRESPONDING  TO  FOUR  DIFFERENT  CONDITIONS  REGARDING  ASPECT 
ANGLES  AND  THE  NUMBER  OF  FREQUENCIES  USED  N. 


W  A&W 


KNOWN  ASP,  N  =  2 


KNOWN  ASP,  N  =  < 


KNOWN  ASP,  N  =  8 


21 

33 

19 

/ 

/ 

25 

17 

23 

17 

/ 

/ 

17 

11 

1 

15 

13 

18 

13 

5 
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Figure  5.12.  Comparison  of  error  probability  between  the  cases 
where  aspect  angles  are  either  known  or  unknown. 
Only  amplitude  is  used.  N  =  8  frequencies. 

A  90%  confidence  interval  for  a  classification 
error  of  30%  is  approximately  ±  5%  . 
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Figure  5.14.  Error  probabilities  for  "known  aspect  angles"  and 
"unknown  aspect  angles".  Both  amplitude  and  phase 
are  used.  Number  of  frequencies  =  8. 

A  90%  confidence  interval  for  a  classification 
error  of  30%  is  approximately  ±  5%  . 


The  table  also  shows  that  increasing  the  number  of  frequencies 
will  significantly  decrease  the  required  S /tf.  In  particular  there 
is  a  decrease  of  20  dB  for  the  'TIME'  algorithm  if  N  is  increased 
from  2  to  8.  It  can  also  be  seen  from  the  table  that  knowledge  of  the 
target  aspect  angle  can  reduce  the  required  S/N  by  about  2  dB. 

E.  CLASSIFICATION  OF  AIRCRAFT 

The  type  of  experiments  discussed  in  section  D  were  applied  to  an 
aircraft  data  base  consisting  of  12  classes,  as  shown  in  Table  5.8  The 
aspect  angles  of  the  aircraft  were  chosen  to  be  0,  90,  and  180  degrees 
so  that  classes  in  the  data  base  are  independent  of  each  other.  The 
aircraft  have  dimensions  of  the  order  of  50  meters,  which,  from  Equation 
(5.27),  implies  that  the  frequency-sampling  interval  for  the  radar 
returns  should  satisfy 

Af  <  3  x  108  =  l  MHz  •  (5.32) 

6  x  50 

Based  on  the  criterion  in  Equation  (5.32)  we  choose  the  frequency¬ 
sampling  interval  Af  as  1  MHz.  The  frequency  ranges  were  chosen  as 
7  MHZ  <  f  <  8  MHZ  for  N  =  2  frequencies,  7  MHZ  <  f  <  10  MHZ  for  N  =  4, 
and  7  MHZ  <  f  <  14  MHZ  for  N  =  8.  The  corresponding  range  of  wavelengths 
for  7  MHZ  <  f  <  14  MHZ  is  21  m  <  x  <  42  m,  and  the  ratio  of  the  aircraft 
dimensions  versus  the  wavelength  is  of  the  order  of 


TABLE  5.8 


A  LIST  OF  AIRCRAFT  USED  IN  THE  CLASSIFICATION  EXPERIMENT 


CLASS 


TARGET 


ASPECT  ANGLE 


Elevation  angle  =  0 


Polarization:  Vertical 


which  means  the  radar  returns  are  in  the  lower  resonance  region.  The 
signal  power  of  the  radar  returns  may  vary  over  a  wide  range.  For  the 
data  base,  the  average  signal  power  was  estimated  to  be  approximately 
equal  to  28  dB-square-meters. 

In  this  case  the  number  of  classes  is  12,  therefore  the  90  percent 
confidence  interval  for  the  probability  of  mi sclassifi cation  Pe  is 
slightly  different  from  that  of  the  ship  classification  where  there  were 
15  classes.  From  Equation  (5.25),  with  M  =  12,  n  =  15,  a  =  0.1,  £ a  = 
1.65,  ?0  ~  0.3,  a  90  percent  confidence  interval  for  Pe  is  given  by 


=  Pe  ±  5.6  percent 

Figures  5.16-5.19  are  plots  of  error  probabilities  for  various 
classification  algorithms  with  the  number  of  frequencies  N  as  a 
parameter.  In  all  four  figures  we  assume  that  target  aspect  angles  are 
known  a  priori,  and  the  number  of  frequencies  N  is  set  to  2,  4,  and  8. 
It  can  be  seen  that  the  error  probabilities  can  be  improved 
significantly  by  increasing  the  number  of  frequencies.  This  is 
especially  true  for  the  time-domain  algorithms. 


Error  probabilities  for  N=2,4,8  frequencies.  Only 
amplitude  data  are  used.  The  aspect  angles  are  assumed 
known.  A  90%  confidence  interval  for  30%  classification 
error  is  approximately  ±  5.6%  . 


Error  probabilities  for  N=2,4,8  frequencies.  Only 
phase  data  are  used.  The  aspect  angles  are  assumed 
known.  A  90%  confidence  interval  for  30%  classification 
error  is  approximately  ±  5.6%  . 
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Figure  5.19.  Error  probabilities  for  N=2,4,8  frequencies.  The 
time-domain  al gv -i thm  is  used.  The  aspect  angles 
are  assumed  known.  A  90%  confidence  interval  for 
30%  classification  error  is  approximately  ±  5.6%  . 
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Figures  5.20-5.22  show  the  effects  of  employing  different  types  of 
classification  algorithms.  In  all  three  figures  we  assume  that  target 
aspect  angles  are  known.  The  number  of  frequencies  N  used  in  each 
figures  are  N  =  8  in  Figure  5.20,  N  =  4  in  Figure  5.21,  and  N  =  2  in 
Figure  5.22.  It  can  be  seen  that  the  performance  of  the  'TIME' 
algorithms  becomes  better  than  the  others  as  the  number  of  frequencies  N 
is  increased.  As  is  the  case  in  ship  classifications,  the  error  in  the 
phase  feature,  W,  dominates  the  N.N.  decision. 

Figure  5.23  shows  the  error  probabilities  when  relative  amplitudes, 
X,  are  employed.  As  is  the  case  for  the  W's,  the  error  in  the  relative 
amplitudes  dominates  the  decisions. 

Figures  5.24-5.27  are  comparisons  of  error  probabilities  between 
the  cases  where  target  aspect  angles  are  either  known  or  unknown. 

Figures  5.24-5.27  show,  respectively,  the  results  of  using  amplitude 
returns  only,  phase  returns  only,  both  amplitude  and  phase  returns,  and 
the  time-domain  algorithms.  Note  if  the  aspect  angle  of  the  unknown 
target  is  known  a  priori,  the  maximum  error  is  4/5  =  80  percent,  and  if 
the  aspect  angle  is  unknown  the  maximum  error  is  11/12  =  92  percent.  It 
is  seen  that  knowledge  of  the  target  aspect  angle  will  improve  the 
results. 

Based  on  the  above  figures  we  can  derive  a  list  of  minimum  signal- 
to-noise  power  ratios  (S/N)  required  for  less  than  10  percent  error 
probabilities  as  shown  in  Table  5.9.  This  gives  us  a  quantitative 
indication  of  the  performance  of  each  algorithm.  From  the  table  we  note 
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that,  for  N  =  2,  the  minimum  S/n  required  for  each  algorithm  has  the 
following  order. 

'A'  <  ' A&W 1  <  'TIME'  <  'W' 

For  N  =  4,  the  order  is  given  as  follows. 

'A'  <  'TIME'  <  'A&W  <  '»' 

And  for  N  =  8  the  order  is  given  as  follows. 

'TIME'  <  'A*  <  'A&W  <  'W*  -  TO'  -  'V 

The  table  also  shows  that  increasing  the  number  of  frequencies  will 
significantly  decrease  the  required  S/a/  .  It  can  also  be  seen  from  the 
table  that  knowledge  of  target  aspect  angles  can  reduce  the  required  S/N 
by  about  4  dB. 
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frequencies  =  8,  and  the  aspect  angles  are  known 
A  90%  confidence  interval  for  30%  classification 
error  is  approximately  ±  5.6%  . 


Error  probabilities  for  four  algorithms.  Number  of 
frequencies  =  8,  and  the  aspect  angles  are  known. 

A  90%  confidence  interval  for  30%  classification  error 
is  approximately  ±  5.6%  . 
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.25.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  Only 
phase  data  are  used,  and  number  of  frequencies  =  8. 

A  90%  confidence  interval  for  30%  classification  error 
is  approximately  ±  5.6%  . 
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Figure  5.27.  Comparison  of  error  probabilities  between  the  cases 
where  aspect  angles  are  known  or  unknown.  The  time- 
domain  algorithm  are  used,  and  number  of  frequencies  = 
A  90%  confidence  interval  for  30%  classification  error 
is  approximately  ±  5.6%  . 


CHAPTER  VI 


CONCLUSIONS  AND  RECOMMENDATIONS 

H.F.  radars  can  detect  targets  at  over-the-horizon  distances,  but 
their  performance  is  dominated  by  ionospheric  conditions  and  sea  states. 
As  a  result,  H.F.  radar  returns  are  usually  accompanied  by  a  large 
amount  of  variability,  and  must  be  extensively  processed.  One  reliable 
processing  technique  is  to  calibrate  the  radar  return  with  respect  to  a 
reference.  For  the  calibration  of  amplitude  returns,  we  can  choose  sea 
scatter  as  a  reference  whose  radar  cross  section  can  either  be 
predicted  or  measured  by  other  means.  We  derived  an  equation  which 
approximately  described  the  relation  between  the  coherent  observation 
time  and  the  post-processing  signal-to-noise  power  ratio  (SNR).  This 
equation  showed  that  we  can  obtain  an  SNR  of  10  dB  for  a  coherent  obser¬ 
vation  time  of  16  minutes.  Higher  SNR  can  be  obtained  for  longer 
coherent  observation  time.  The  measurement  of  the  phase  of  the  target 
returns  requires  an  elaborate  processing  scheme  because  the  phase  path 
length  may  vary  over  a  wide  range,  resulting  in  ramp- type  phase  errors 
as  a  function  of  frequency. 

Because  of  the  errors,  the  phase  of  the  measured  target  returns  in 
most  cases  cannot  be  directly  chosen  as  a  set  of  features  for  classifi¬ 
cation.  Instead  we  derived  a  variable  W  from  the  phase,  and  showed  that 
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the  W's  were  Independent  of  wave  path  length  and  therefore,  did  not 
require  a  phase  reference.  The  amplitude  of  the  RCS  and  a  set  of  nor¬ 
malized  W's  formed  a  feature  set  for  the  frequency-domain  Nearest- 
Neighbor  classification.  The  Nearest-Neighbor  decision  rule  assigned  a 
target  to  a  class  whose  distance  from  the  target  was  the  smallest  among 
all  classes  in  the  data  base.  Time-domain  data  were  obtained  by  taking 
the  Discrete  Inverse  Fourier  Transform  (DIFT)  of  the  multi -frequency 
measurements  of  the  target  returns.  The  DIFT  converted  the  ramp-type 
phase  errors  to  shifts  of  waveforms  along  the  time  axis.  By  choosing 
the  maximum  value  of  a  cross-correlation  function  of  two  waveforms, 
we  essentially  eliminated  the  difference  in  time-shifts  between  the 
waveforms.  The  time-domain  decision  rule  assigned  a  target  to  a  class 
whose  maximum  cross-correlation  with  the  target  was  the  largest  among 
all  classes  in  the  data  base. 

A  set  of  radar  backscatter  measurements  of  scale  models  of  ships 
and  aircraft  was  used  to  generate  a  data  base  for  testing  the  target 
classification  algorithms.  The  measurements  were  conducted  in  an 
anechoic  chamber  at  certain  scaled  frequency  ranges,  e.g.,  1-2GHZ, 
2-4GHZ,  4-8GHZ,  or  8-12GHZ.  Connecting  some  of  these  frequency  bands, 
and  removing  possible  phase  offsets,  we  obtained  a  number  of  data  files 
of  wider  frequency  ranges.  The  resulting  data  files  were  convolved  with 
a  Hamming  window  to  remove  clutter  present  in  the  waveforms,  and,  at 
the  same  time,  estimate  the  values  of  the  data  at  selected  frequencies. 
The  Hamming  window  was  equivalent  to  a  time-domain  low-pass  filter  with 


a  properly-designed  time-width,  and  its  center  coincided  with  the 
target's  center  so  that  it  effectively  removed  clutter  returns  located 
outside  the  target  zone.  The  filtered  results  were  scaled  according  to 
scale  factors  for  each  target  and  a  full-scale  data  base  was  generated. 

The  testings  of  the  classification  algorithms  made  use  of  the  com¬ 
puted  full-scale  data  base,  and  were  conducted  using  FORTRAN  programs  on 
a  DEC  VAX  computer.  By  adding  Gaussian  noise  to  each  class  in  the  data 
base,  we  created  a  set  of  test  target  data  upon  which  the  classification 
algorithms  were  applied.  The  classification  experiments  were  repeated 
several  times,  and,  based  on  the  results,  a  Maximum  Likelihood 
Estimator  (MLE )  was  applied  to  estimate  the  probabilities  of  misclassi- 
fi cation.  We  constructed  two  data  bases  for  testings  of  both  ship  and 
aircraft  classification  techniques.  The  former  consisted  of  15  indepen¬ 
dent  classes,  corresponding  to  5  ships  at  0,  90,  and  180  degrees  aspect 
angles,  and  the  latter  consisted  of  12  classes  corresponding  to  4  air¬ 
craft  at  0,  90  and  180  degrees  aspect  angles.  The  frequency  ranges  were 
chosen  as  1.7  MHz  <  f  <  4.5  MHz  for  the  ships,  and  7  MHz  <  f  <  14  MHz 
aircraft  so  that  the  ratios  of  the  target  dimensions  to  the  wavelengths 
were  of  the  order  of  1,  which  implied  that  the  radar  frequencies  were 
located  in  the  lower  Resonance  region.  Based  on  the  Shannon's  sampling 
theorem,  the  frequency  increments  were  selected  as  1  MHz  for  the  ships 
and  0.5  MHZ  for  the  aircraft.  The  aspect  angles  were  assumed  known  for 
most  experiments  in  order  to  obtain  lower  errors.  We  compared  the  mis- 


classification  percentages  as  a  function  of  post-processing  signal-to- 
noise  power  ratio  (SNR)  with  respect  to  the  following  cases. 

1.  Use  of  2,  4,  or  8  frequencies. 

2.  Use  of  different  classification  algorithms  including  the 
frequency-domain  and  the  time-domain  algorithms. 

3.  Target  aspect  angles  are  known  or  unknown  a  priori. 

From  the  results  of  case  1  we  found  that  increasing  the  number  of 
frequencies  would  significantly  decrease  the  error  probabil ities, 
especially  for  the  time-domain  curves.  From  the  results  of  case  2,  we 
found  that  for,  a  small  number  of  frequencies  the  '  amplitude  only  1 
algorithm  had  the  lowest  error  probabilities,  while,  for  high  number  of 
frequencies,  the  'time-domain'  algorithm  had  the  lowest  error  probabil¬ 
ities.  The  results  showed  that  for  an  SNR  of  greater  than  10  dB, 
which  means  that  the  noise  power  is  about  10%  of  the  signal  power,  it 
would  be  possible  to  obtain  a  mi sclassifi cation  error  of  less  than  10%. 
The  'VI'  and  the  'relative  amplitude'  algorithms  had  the  highest 
mi sclassifi cation  errors.  It  was  noted  that,  when  two  incompatible 
features  such  as  'amplitude  and  W  or  'relative  amplitude  and  W  were 
to  be  used  together,  they  should  be  properly  weighted  in  order  to 
get  better  results.  From  the  results  of  case  3  we  found  that  a  priori 
knowledge  of  the  aspect  angle  would  decrease  the  error  probabilities, 
and  not  by  as  much  as  an  increase  in  the  number  of  frequencies. 


i 

This  study  shows  that  classifications  of  targets  at  over-the- 
1  horizon  distances  by  processing  H.F.  radar  returns  is  feasible.  For  a 

reasonable  coherent  observation  time  such  as  16  minutes  or  longer,  we 
can  obtain  an  SNR  of  about  10  dB  or  higher,  which,  according  to  our 
experimental  results,  supports  the  possibility  that  we  can  classify  the 
targets  with  less  than  10  percent  errors.  We  conclude  that 
classification  of  targets  using  HF  radars  is  feasible. 

It  is  recommended  that  future  work  be  focused  on  five  areas.  First 
of  all,  in  order  to  minimize  the  probability  of  misclassification  it 
is  recommended  that  optimum  normalization  constants  for  both  the  phase 
features  W  and  the  relative-amplitude  features  X  be  found.  An  optimum 
range  of  frequencies  which  minimizes  the  probabilities  of  misclassi¬ 
fication  should  also  be  found.  Secondly,  a  larger  data  base  should 
be  established.  Even  though  the  above  data  base  served  to  indicate 
classifiability,  it  was  not  large  enough  to  draw  adequate  statistical 
conclusions  on  the  expected  probabilities  of  misclassification. 

Thirdly,  uncataloged  targets  should  be  considered.  Fourthly,  since 
previous  efforts  to  classify  targets  were  based  upon  a  set  of  data 
obtained  from  measurements  of  scaled  model  targets.  It  is  recommended 
that  data  obtained  from  full-scale  OTH  radar  systems  be  used  to  test 
the  classification  algorithms  developed  here,  and  to  provide  insight 
into  new  approaches  for  classification.  And  lastly,  methods  for 
measuring  more  accurately  the  phase  and  amplitude  of  the  target  returns 
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APPENDIX  A 


This  appendix  contains  plots  of  amplitude  and  phase  of  the  pro¬ 
cessed  radar  returns  of  five  ships  and  four  aircrafts.  The  aspect 
angles  are  0,  90,  and  180  degrees. 
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Figure  A.l.  Ship  A  at  0  degrees 
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Figure  A. 22.  Aircraft  H  at  0  degrees 


APPENDIX  B 


The  computer  programs  developed  for  this  research  are  listed  in 
this  appendix  with  comments.  The  program  names  are  DOC. FOR,  133. FOR, 


flD-fll 62  449 


AUTOMATIC  TARGET  CLASSIFICATION  USING  HF  MULT I FREQUENCY 
RADARS (U)  OHIO  STATE  UNIV  COLUMBUS  ELECTROSCIENCE  LAB 
J  CHEN  OCT  83  ESL-714190-2  N00014-82-K-0037 


3/3 


UNCLASSIFIED 


F/G  17/9 


NL 


PROGRAM  NAME  :  DOC. FOR 


i 

I 


C 
C 

c 

C  THIS  PROGRAM  CALIBRATES  MEASURED  TARGET  RETURNS 

C 

C 

REAL*4  AM1201 ) ,PH(201 ) 

INTEGER *2  TFILEI6) .BFILE16) ,CFILE(6) ,EFILE(6) , RFILE16) 
INTEGER*2  BUFFd024) 

CHARACTER *12  RSFILE 

COMPLEX*  B  TA ( 201J , B A ( 201 ) , CAS ( 201 1 , EXS ( 201 ) , RS { 201 ) 
DATA  N/201/ 

C 

C  SET  UP  AN  OUTPUT  BUFFER 

C 

COMMON  /BLK1/BUFF 
C 

C  FILE  NAMES  ARE  STORED  IN  INTEGER  ARRAYS  AND 

C  ARE  INITIALLY  SET  TO  0. 

C 

288  DO  147  1=1,6 

TFILE(I 1=0 
BFILE(I]=0 
CFILE(I )=0 
EFILEd  1=0 
147  RFILEd  1=0 

10  F0RMAT11X,/,'  INPUT  THE  TARGET' 1 

C 

C  READ  A  TARGET  FILE  FROM  A  STORAGE  UNIT 

C  l  IN  MOST  CASES  FROM  A  FLOPPY  DISK  3 

C 

TYPE  10 

ACCEPT  120.TFILE 
120  FORMAT (6A2) 

CALL  T DATA! RSFILE, TFILE.TA.O] 

C 

C  READ  A  BACKGROUND  FILE 

C 

TYPE  20 

20  FORMATdX,/,'  INPUT  THE  BACKGROUND' ) 

ACCEPT  120.BFILE 
CALL  T DATA ( RSFILE, BFILE, BA, 01 
C 

C  READ  A  CALIB RATION-SPH ERE  FILE 

C 

TYPE  40 

40  FORMATdX,/,'  INPUT  THE  CALIBRATE  SPHERE' 1 
ACCEPT  120.CFILE 
CALL  T DATA ( RSFILE, CFILE, CAS, 0) 

C 

C  READ  A  THEORETICAL  CALIBRATION-SFHERE  FILE 

C 

TYPE  50 

50  FORMATdX,/,'  INPUT  THE  EXACT  SPHERE') 

ACCEPT  120.EFILE 


180 


CALL  TQATAfRSFILE.EFILE.EXS.O) 
CALIBRATE  THE  TARGET  FILE 


IF  THE  CALIBRATED  RESULTS  ARE  IDENTICAL 
TO  ZBU)  SET  THEM  TO  1.E-15 

DO  260  1=1 ,N 

BS  tl )  =  tTA(I ) — B A ( I ) ) • EXS ( I ) / { CAS [ I ) — B A ( I ) ) 

IF (AIMAG1 RS ( I ) ) . EQ. O.AND. REAL(RS(I ) ) .EQ.O. ) 

RS(I )=CMPLX(1 .E-15,1 . E— 1 5 ) 

CONTINUE 

STORE  THE  FILE  NAMES  IN  THE  OUTPUT  BUFFER 

DO  478  1=1 ,6 
BUFF (24+1 )=TFILE{I) 

BUFF  ( 24+6+1 )  =8FILE  II) 

BUFF  1 24+1 2+1 1  =CFILE  (I) 

BUFF  ( 24+1  B+I )  =EFILE  (I) 

CONTINUE 

STORE  THE  CALIBRATED  TARGET  FILE  ON  A  STORAGE  UNIT 

CALL  TDATA(RSFILE,RFILEtRS,1 ) 

TYPE  DO  YOU  WANT  TO  PROCEED?  Y=1  M=0' 

ACCEPT  *, IYN 
IF1IYN.ED.1 ) GO  TO  2BB 
CALL  EXIT 
END 


SUBROUTINE  TDATAIRFILE, INFILE, CA.IRW) 

THIS  SUBROUTINE  HANDLES  DATA  I/O 
SET  IRW  =  0  FOR  INPUT 
SET  IRW  =  1  FOR  OUTPUT 

COMPLEX’S  CA (201 ) 

INTEGER *2  BUFF(1024) ,10(48) ,IC0DE(4] 

REAL *4  AM (201 ) ,  PH  (201 ) 

INTEGER*2  INFILE(6) 

CHARACTER *12  RFILE 

EQUIVALENCE  (10(1) .BUFF ( 1 ) ] , ( ICODE ( 1 ) .BUFF ( 49) ) , 
(AM(1),BUFF(53)),(PH(1), BUFF (455)) 

COMMON  /BLK1/BUFF 
IFtlRW.BQ.I  )G0  TO  10 

THIS  SECTION  READS  A  FILE 

THE  FILES  TO  BE  READ  WERE  ORIGINALLY 
WRITTEN  ON  DISKS  IN  R-T  FORMATS  ON  A  PDP-11/03 
COMPUTER,  WHICH  MONITORED  THE  MEASUREMENTS. 

THE  FILE  FORMATS  WERE  TRANSFORMED 


TO  VAX  FORMATS,  AND  CAN  BE  READ  BY 
A  VAX  COMPUTERS. 

THE  TRANSFOMATION  CONVERTS  THE  DATA 
INTO  ASCII  CODES,  AND  STORED  THE 
THE  CODES  IN  A  256-INTEGER  BUFFER, 

WHERE  EACH  INTEGER  CAN  HOLD  TWO  ASCII  BYTES. 


OPEN  A  FILE  WHOSE  FILE  NAME  IS  STORED  IN 
THE  ARRAY  INFILE 

OPEN ( UNIT=1 9 , NAME=I NFILE, TYPE= ' OLD ' ) 

DO  20  1=1,1024,256 
J=I+255 

READ11 9,30]  (BUFF(K)  ,K=I,  J) 

FORMAT! 256 A2) 

THE  UNITS  OF  THE  AMPLITUDE  RETURNS  STORED  IN  THE 
ORIGINAL  FILES  ARE  OB-SQUARE-CM,  ANO  THE  UNITS  OF 
THE  PHASE  RETURNS  ARE  RADIANS. 

THE  FOLLOWING  CONVERTS  THE  DATA  INTO  COMPLEX  FORM 

DO  70  1=1,201 

AM(I  )=10**(AM(I  )/20.  ] 

CA(I)=CMPLX(AM(I )*COS(PH (I } ) ,AM(I )*SIN(PH(I) ) J 
GO  TO  BO 

THIS  SECTION  WRITES  A  CALIBRATED  TARGET 
RETURNS  ON  A  STORAGE  UNIT 

CONVERTS  THE  COMPLEX  DATA  INTO 
AMPLITUDE-PHASE  FORMAT 

DO  210  1=1,201 

AM(I ]=20.*AL0G10(CABS(CA(I ) ] ) 

PH ( I ) =ATAN2 ( AI MAG ( CA ( I )J , REAL ( CA (I) ) ] 

SPECIFY  A  NAME  FOR  THE  OUTPUT  FILE 

TYPE  92 

ACCEPT  94.RFILE 

WRITE  A  4 B- CHARACTER  SENTENCE  EXPLAINING  THE  OUTPUT  FILE 
TYPE  112 

F0RMAT11X,/,'  INPUT  THE  ID  FOR  THE  RESULTS') 

ACCEPT  128,  ( I D ( I )  ,1=1 ,24} 

FORMAT (24A2) 

FORMAT  (IX,/, '  INFUT  A  FILE  NAME  FOR  THE  RESULTS') 

FORMAT (AI 2) 

OPEN  AN  OUTPUT  FILE 

OPEN ( UNIT=1  9, NAME=RFILE,TYPE= '  NEW ' ,  RECORDS  12  E=51 2 ) 


uuuuuuu 


PROGRAM  NAME  :  B3.P0R 


THIS  PROGRAM  MAKES  USE  OF  THE  2-4GH2,  4-8GHZ,  AND  8-12GHZ 
CALIBRATED  DATA  FILES  GENERATED  BY  THE  PROGRAM 
DOC. FOR,  AND  GENERATE  A  FULL-SCALE  DATA  FILE. 


REAL*4  FA  ( 801 ) ,  SA  ( 801 ) ,  A24(  201 ) ,  A48(  201 ) ,  A  82  ( 201 ) 

1  .P24C201 ) ,P48(201 ) ,P82(2D1 ] ,FA24(201 ) ,FA4B(201 ) ,FAB2(201 ] 
INTEGER*4  10(12) 

INTEGER*2  FI 2(6) ,F24(6] ,F48(6) ,F82(E) ,RFILE(6) ,LINE2(24) 
COMPLEX* 8  TA(201 ) , CA( 801 ) ,R( 801 ) 

INPUT  THE  SCALE  FACTOR  OF  A  TARGET,  A  DESIRED 
FREQUENCY  INCREMENT  FOR  THE  OUTPUT  FILE,  AND 
THE  WIDTH  OF  THE  HAMMING  WINDOW. 

TYPE  *, 'INPUT  THE  SCALE  FACTOR,  FREQ  INC  AND  WINDOW  SIZE  (MHZ) 
ACCEPT  *,SF,FINC,WD 

READ  THE  NAME  OF  A  2-4GHZ  CALIBRATED  DATA  FILE 
TYPE  20 

FORMAT (IX,/, '  INPUT  THE  2-4G  FILE  ') 

ACCEPT  120.F24 
FORMAT (6A2) 

BEAD  THE  NAME  OF  A  4-8GHZ  CALIBRATED  RETURNS 
TYPE  40 

FORMAT (IX,/ , 1  INPUT  THE  4-BG  FILE') 

ACCEPT  120.F48 

READ  THE  NAME  OF  A  8-12GHZ  CALIBRATED  RETURNS 
TYPE  50 

FORMAT (IX,/, '  INPUT  THE  B-12G  FILE') 

ACCEPT  120.F82 

SPECIFY  A  NAME  FOR  THE  OUTPUT  FILE 

TYPE  *,'  INPUT  A  FILE  NAME  FOR  THE  RESULT' 

INPUT  A  FILE-HEADER  WHICH  EXPLAINS  THE  OUTPUT 
FILE 


ACCEPT  120.RFILE 

TYPE  *,'  INPUT  A  HEADER  (48  CHARACTERS)' 
ACCEPT  732,  (ID(IZ)  ,IZ=1 ,12) 

FORMAT (12A4) 

INITIALIZE  THE  FILE-NAME  ARRAYS 


F24(6)=0 

F24(5)=0 

F48(6)=0 

F48I5)=0 

F8216)=0 

F  82 [51 =0 

RFILE(6)=0 

RFILE(5)=0 

CALCULATE  FULL-SCALE  FREQUENCIES 
ACCORDING  TO  THE  SCALED  FACTOR 

DO  46  1=1,201 

FA241I )=(2.+(I-1 )/100.)/SF*1000. 
FA48(I)=(4.+(I-1)/50.)/SF*1000. 

FA82(I )=t  8.+II-1 )/50. )/SF*1 000. 

READ  THE  2-4GHZ,  4-8GHZ,  AND  8-12  GHZ  FILES 

CALL  T0ATA(F24,A24,P24,P24AIP24Z) 

CALL  TDATAIF48,A48,P48,P4BA,P4BZ] 

CALL  TDATA(F82,AB2,PB2,P82A,P82Z) 

CALCULATE  PHASE  OFFSETS  AT 

THE  JUNCTIONS  OF  THE  FILES,  NAMELY,  AT  4  GHZ, 

AND  B  GHZ. 

DPL=P24Z-P4  8A 
DPR=P  82A-P4  BZ 

CORRECT  THE  PHASE  OFFSETS  IN  THE 
2-4GHZ,  AND  8-12  GHZ  FILES, 

AND  USE  THE  4-BGHZ  FILE  AS  REFERENCE 

CALL  PHCR  ( DPL ,  FA24, 201 ,  P24] 

CALL  PHCR ( DPR , FA82 , 1 ,PB2) 

CONNECT  THE  CORRECTED  FILES  (  2-4GHZ,  4-BGHZ,  AND  8-12  GHZ  ) 
AND  CREATE  A  FULL-SCALE  COMPLEX  TARGET 
RETURNS  ACCORDING  TO  THE  SCALE  FACTOR 

DO  3540  1=1,201 

CA (I) =CHPLX ( SF*A24 1 1 1  *COS l P24 ( 1 1 1 , SF* A24CI) *S I N ( P24( I ) ) 1 
CA[I+200)=CMPLX(SF*A48(I )*COS  (P48(I) )  ,SF*A4B(I )«SIN(P4B(I ) ) ) 
CA ( I +400 ] =CMPLX ( SF*A82 l I ) *COS ( P 82 ( I JJ , SF* A82 l I ) *SI N { P 82 ( I ) ) ) 
SA(I 1=2. 

SA(I+200. )=4. 

SAII+400. 1=4. 

FA(I)=FA241I ) 

FAII+200)=FA4B(I ) 

FA(I+400)=FA82(I) 

STORE  THE  INRJT-FILE  NAMES  IN  A  BUFFER 


987  LINE2d  )  =  '  ' 

DO  999  1=1,4 
LINE2(I+6!=F24(I) 

LINE2(I+12]=F4B(I ) 

LINE2II+1  8)=FB2d) 

999  CONTINUE 

C 

C  SMOOTH  THE  FULL-SCALE  DATA  FILE  BY  CONVOLVING 

C  THE  FILE  WITH  A  HAMMING  WINDOW 

C 

C  CALCULATE  PARAMETERS  FOR  THE  CONVOLUTIONS,  WHERE, 

C  FLOW  IS  THE  LOWEST  FULL-SCALE  FREQUENCY 

C  FHIGH  IS  THE  HIGHEST  FULL-SCALE  FREQUENCY 

C  FINC  IS  THE  FULL-SCALE  FREQUENCY  INCREMENT 

C 

FLOW=2./SF*1000. 

FHIGH=1 2 ./SF*1 OOO , 

FLOW=I NT { FL0W*1 0 . ) /1 0 .  + . 1 
FHIGH=INT(FHIGH*10.)/10. 

1=1 

K=0 

ICNT=601 

TYPE  * , ’ FA ( 1 )  =  * , FA C 1 ) » 1  FA{ICNT]  =  ' .FAtICNT) 

TYPE  *,’SF=',SF,'  ICNT=' , ICNT, 1  FLOW=' .FLOW, '  FHIGH=’ , FHIGH 
C 

C  PICK  UP  THE  DESIRED  FREQUENCY  AND  POSITION  THE 

C  CENTER  OF  THE  HAMMING  WINDOW  AT  THAT  FREQUENCY. 

C  ESTIMATE  THE  DATA  VALUE  OF  THE  DESIRED  FREQUENCY  BY 

C  TAKING  THE  WEIGHTED  AVERAGE  OF  THE  NEIGHBORING  DATA  POINTS 

C  COVERED  BY  THE  HAMMING  WINDOW,  WITH  THE  WEIGHTINGS  DETERMINED 

C  THE  HAMMING  WINDOW. 

C 

DO  7720  RF=FLOW, FHIGH, FINC 
7740  IFtRF.GE.FAd] .AND.RF.LE.FA(I+1 ) )GQ  TO  7730 

1=1+1 

IFd.GT.ICNTIGO  TO  4912 
GO  TO  7740 
7730  K=K+1 

CALL  INTER(R,CA,FA,SA, ICNT, NS, I , K, WD, RF) 

7720  CONTINUE 

C 

C  OUTPUT  THE  SMOOTHED  FULL-SCALE  DATA  FILES 

C 

4912  CALL  P0ATA(RFILE,R,K,FL0W,FINC,SF,I0,LINE2) 

C 

C 

TYPE  NPT=',K, *  FLOW=',FLOW, '  FINC=',FINC 
TYPE  *, 1  DO  YOU  WANT  TO  PROCEED  ?  Y=1,N=0' 

ACCEPT  *,IYN 
IFdYN.EQ.DGO  TO  3542 
CALL  EXIT 
END 


SUBROUTINE  T DATA (INFILE, AMT, PHT,PHA,PHZ! 


I 

i 


> 


c 

C  THIS  SUBROUTINE  READS  A  DATA  FILE 

C 

INTEGER*2  BUFF(1024)  ,INFILE(6] 

REAL*4  AM ( 201 ) , PH ( 201 ) , AMT( 201 ) , PHT 1 201 ) 

EQUIVALENCE  (AM(1) , BUFFI 53 ) ),(PH (1) ,BUFF ( 455) ) 
0PEN(UNIT=1 9,NAME=INFILE,TYPE='0L0' ) 

PI=4.*ATAN(1.) 

TPI=PI*2. 

00  20  1=1,1024,256 
J=I+255 

20  READ(19,30)  [BUFF C K ]  ,K=I,  J) 

30  FORMAT! 256 A2) 

DO  10  1=1,201 
I F  CAM  II) .GT.40. )AM[I )=40 . 

AMT(I)=10**(AM(I)/20.) 

PHT(I )=PH (I) *1  BO. /PI 
10  CONTINUE 

C 

C  ELIMINATE  BRANCH  CUTS  IN  THE  PHASE  RETURNS 

C 

CALL  LPV(PHT,1 B5. ,30. ,1 90. ,-1 80. , 201 ) 

C 

DO  55  1=1,201 

55  PHT  1 1 )  =PHT  ( I )  *PI/1 80 . 

CL0SE{UNIT=1 9,DISP='SAVE' ) 

C 

C  ESTIMATE  THE  PHASE  VALUES  AT  THE  ENOS  OF  THE  FILE 

C 

CALL  EPH[PHA,PHT,2,7,1) 

CALL  EPHIPHZ, PHT, 200, 195,-1) 

C 

RETURN 

END 

C 

SUBROUTINE  POATACRFILE.CA, NPT, FLOW, FINC, SF, ID, LINE2) 
C 

C  THIS  SUBROUTINE  WRITES  A  FILE  ON  A  STORAGE  UNIT 

C 

INTEGER*2  RFILE16) ,LINE2C24) 

INTEGERS  ID(12] 

REAL*4  AMI  801), PH  (801) 

COMPLEX* 8  CA1801) 

10  DO  210  1=1,  NPT 

IF(CABS(CA(I ) ) .EQ.O. )CA(I )=(1 . E-15, 1 . E-1 5) 

AM(I )=20 . *AL0G1 OICABStCAtl ) ) ) 

210  PH(I  )=ATAN2(AIMAG(CA  tl)  J ,  REAL  CCA  [I ) ) ) 

OPEN (UNIT=1 9, NAME=RFILE , TYPE= ' NEW 1 , RECORDSIZ  E=51 2 ) 
WRITEC1 9,30)  (10(1),  1=1, 12),  LINE2,  NPT,  FLOW, FINC, SF 
30  FORMAT (12A4,24A2,4A4) 

WRITE (1 9,40) (AM(I ) ,1=1 ,NPT) , (PH  (I ) ,1=1 , NPT) 

40  FORMAT  (12BA4) 

80  CLOSE (UNIT=1 9, OISP=' SAVE' ) 

RETURN 

END 
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THIS  SUBROUTINE  CALCULATES  THE  WEIGHTED  AVERAGE 
OF  THE  DATA  POINTS,  USING  A  HAMMING  WINDOW  TO  DECIDE 
THE  WEIGHTINGS 

SUBROUTINE  INTER (R, CA, FA, SA, ICNT, NS, I, K, WD, RF) 
COMPLEX’S  RtaOD.CAlBOl) 

REAL*4  FA(801),SA(B01] 

WEI=0 

IS=I 

XTMP=0. 

YTMP=0. 

RFL=RF-WD/2. 

RFH=RF+WD/2. 

IF (IS. GT. ICNT) GO  TO  10 
IF(FA(IS) .GT.RFHJGO  TO  10 
T=FA(IS)-flF 

H AMM= . 54+ . 46 ’COS ( 3 . 1 41 5 926 *T/WD 1 
XTMP=XTMP+HAMM*SA ( IS ) ’REAL ( CA ( IS ) ) 

YTMP= YTMP+HAMM’SA l IS ) ’AIMAG ( CA { IS ) ) 

WEI=WEI +HAMM*SA ( IS ) 

IS=IS+1 
GO  lO  20 
IS=I-1 

IFIIS.LT.I )G0  TO  30 
IF(FA(IS)  .LT.RFDGO  TO  30 
T=FA(I5)-RF 

HAMM= .  54+ .  46*C0S  ( 3 . 1 41 5  926  ’T/WO ) 

XTMP=XTMP+HAMM*SA (IS ) ’HEAL ( CA (IS ) ) 

YTMP=  YTMP+HAMM’SA  {IS )  *AIMGt CA  [IS) ) 

WEI=WEI+HAMM*SA (IS ) 

IS=IS-1 
GO  TO  15 

R (K]=CMPLX(XTMP/WEI , YTMP/WEI 1 

RETURN 

END 

THIS  SUBROUTINE  CORRECTS  THE  PHASE  OFFSETS 

SUBROUTINE  PHCR(DPH,FA,IBASE,PH) 

REAL*4  PH (201), FA (201) 

TW0PI=8.*ATAN(1 . ) 

PI=TW0PI/2. 

IFtDPH .LT.PI )G0  TO  20 
DPH=DPH-TWOPI 
GO  TO  30 

IF(DPH.GT.-PI )G0  TO  40 
DPH=OPH+TWOPI 
GO  TO  20 
FB ASE=FA ( IB ASE ) 

DO  10  1=1 ,201 

PH  (I )  =PH  (I ) -(DPH’FA (I  l/FBASE ) 

RETURN 


C  THIS  SUBROUTINE  ESTIMATES  THE  VALUES  OF 

C  PHASE  RETURNS  AT  THE  ENDS  OF  A  DATA  FILE, 

C  USING  A  LINEAR  INTERPOLATION 

C 

SUBROUTINE  EPH(PHM,PH,NS, NF,NI) 

REAL*4  PH (201) 

X=0 

Y=0 

XY=0 

X2=0 

DO  10  I=NS,NF,NI 
X=X+FLOAT(I ) 

Y=Y+PH (I ) 

XY=XY+FLOAT { I ) *PH ( I ) 

10  X2=X2+FL0AT(I )*FLOAT (I ) 

XN=(NF-NS)*NI+1 

DELTA=XN*X2-X*X 

A= (XN*XY-X*Y)/DELTA 

B=(X2*Y-X*XY) /DELTA 

PHM=A*NS+fi 

RETURN 

END 
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nno  non  nooo  nooo  oonnooon 


PROGRAM  NAME  :  FREQ. FOR 


THIS  PROGRAM  CLASSIFIES  TARGETS  IN  THE 
FREQUENCY  DOMAIN. 


COMPLEX  C(1 8,50) , NC(1 8,50 J 

REAL  D(1  8, 1 8) , A(1 8,50) ,P(1 8, 50] ,W(1  8,5D) , AN(1  8,50) , 

1 PN  ( 1  8, 50 ) ,  WN  ( 1  8, 50 ) 

REAL  WL ( 50 ) , PX{30) ,PY(10,30) ,PYT(30) , SAP [1 8) 

TYPE  */  *****  YOU  MUST  ASSIGN  DATA  DIRECTORY  TO  R  *•***' 

TYPE  WANT  TO  PRINT  DISTANCES  IN  CLASSIFICATION  ?  Y=1  N=0' 
ACCEPT  *, IDIS 

INPUT  TWO  ARBITRARY  LARGE  ODD  INTEGERS 

WHICH  WILL  SERVE  AS  SEEDS  FOR  GENERATEING  RANDOM  NUMBERS 

TYPE  *, 'INPUT  2  ARBITRARY  LARGE  ODD  INTEGERS' 

ACCEPT  *,IS3,IS4 

SPECIFY  THE  MINIHIM,  MAXIHIM,  AND  INCREMENT  OF 
ADDITIVE  GAUSSIAN  NOISE  POWERS 

TYPE  *, ’INPUT  MIN,  MAX,  INC  OF  THE  NOISE  POWER  (DB)’ 

ACCEPT  *,IV1 ,IV2,IVT 

SPECIFY  TARGET  ASPECT  ANGLES 

TYPE  ‘/INPUT  MIN, MAX, INC  OF  ASP1 
ACCEPT  *, MINASP, MAXASP, INCASP 

SPECIFY  FREQUENCIES 

TYPE  ‘/INPUT  FMIN, FINC, NO  OF  FREQ’ 

ACCEPT  *, FMIN, FINC, NF 

TYPE  */  DO  YOU  WANT  TO  USE  RELATIVE  AMPLITUDE  ?  Y=1  N=0' 
ACCEPT  *,  IRA 

TYPE  ’/INPUT  NO  OF  EXPERIMENTS' 

ACCEPT  *,  NEX 

TYPE  *,  '  ASP  KNOWN  - >  TYPE  1  ' 

TYPE  *,'  ASP  UNKNOWN  - >  TYPE  O' 

ACCEPT  *,  IAK 

FMAX=FMIN+(NF-1 )*FINC 

PRINT  */  LARGE  INTEGERS  :  \IS3,IS4 

PRINT  */  AMPLITUDE  AND  PHASE' 

IFtIRA. EQ. 1 1PRINT  *, '  RELATIVE  AMPLITUDE  USED' 

IF  (IAK.  HI.  1 )  PRINT  *,'  ASP  ASSUMED  KNOWN' 

IFdAK.  EQ.  0 1PRINT  *,'  ASP  ASSUNEO  UNKNOWN' 

PRINT  */  ADDITIVE  NOISE  ADDED' 

PRINT  */  CLASSIFICATION  IN  FREQUENCY  DOMAIN' 

PRINT  *,'  NO  OF  EXPERIMENTS  \NEX 
PRINT  2  83  9, FMIN,  FMAX,  FINC,  NF 
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2839  FORMAT ( 1  FMXN  =  ',F8.3,'  FMAX  =  ',F8.3,'  FINC  =  1 
1F8.3,/, '  NO  OF  FRED  =  ',13) 

PRINT  2B49,MINASP,MAXASP,INCASP 
2849  FORMAT ( 1  MIN,MAX,ANO  INC  OF  THE  ASPECTS  :  ',315) 

C 

C  READ  DATA  FILES 

C 

CALL  OATAIMINASP, MAXASP, INCASP, FMIN, FMAX, FINC, 
INF, NS, C, A, P, IRA) 

C 

C  COMPUTE  WAVELENGTHS 

C 

6664  DO  99237  1=1  ,NF 

99237  WL(I )=300./ (FMIN+lI-1 )*FINC) 

PRINT  *, 1  NO  OF  TARGETS  '  ,NS 
C 

C  COMPUTE  W(J)  AND  RELATIVE  AMPLITUDES 

C 

PI=4.*ATAN(1 .) 

TWOPI=2.*PI 
NF1  =NF-1 
DO  1434  1=1  ,NS 
DO  1434  J=1  ,NF-1 

IFIIRA.EQ.1)A{I,J)=A(I,J)/A(I,J+1) 

C 

C  ELIMINATE  BRANCH-CUTS 

C 

PC=P(I,J) 

DP=P(I,J)-P(I,J+1] 

I F ( DP . GT . PI ) PC=PC~TWOPI 
I F ( DP . LT . -PI ) PC=PC+TWOPI 
C 

1434  W(I,J]=WL(J)*PC-WL(J+1)*P(I,J+1] 

C 

C  ESTIMATE  VARIANCES  OF  THE  AMPLITUDES 

C 

IF (IRA. HI. 0) CALL  VARN (A, NF, NS, VA, 1 ) 

IF (IRA. EQ. 1 1CALL  VARN (A, NF1 , NS, VA, 1 ) 

C 

C  ESTIMATE  VARIANCES  OF  THE  W’S 

C 

CALL  VARN (W, NF1 , NS, VW, 1 ) 

C 

PRINT  *,'VAR(AMP)=',VA 
PRINT  * , 1  VAR ( W )  = ' , VW 
C 

C  ESTIMATE  SIGNAL  POWERS 

C 

2880  PAVE=0 . 

DO  910  1=1,  NS 
DO  910  J=1  ,  NF 
PAVE=PAVE+A(I,J)*A(I,J) 

910  CONTINUE 

PAVE=PAVE/ ( FLOAT (NF)*FL0AT( NS) ) 

PAVE=1 0 . *AL0G1 0 ( PAVE ) 


PRINT  V  AVE  SIGNAL  P0WB1  (OB)  =  '.PAVE 

DEFINE  THE  VARIANCE  RATIO 

TYPE  *, 'INPUT  R=VAR (Al/VAR (KW) ' 

ACCEPT  *,RR 

SPECIFY  THE  FEATURES  AS  FOLLOWS 

IYW  =  1  - >  USE  THE  W'S 

IYW  =  0  - >  DO  NOT  USE  THE  W'S 

IYA  =  1  - >  USE  THE  AMPLITUDES 

IYA  =  0  - >  DO  NOT  USE  THE  AMPLITUDES 

TYPE  *, 'INPUT  IYW  AND  IYA' 

ACCEPT  *, IYW, IYA 

CALCULATE  THE  NORMALIZATION  CONSTANT  RK 
RK=SQRT(VA/(RR*VW) ) 

PRINT  18220 

FORMAT (IX  ' ***************************************  < ) 

PRINT  *, ' IYW=' ,IYW, ' IYA=' ,IYA 
PRINT  18221  ,RR 

F0RMAT(1X, ' VAR(A]/VAR(KW]=' .E12.5) 

VARY  THE  NOISE  POWER  FROM  IV1  TO  IV2  OB 
AT  AN  INCREMENT  OF  IVT  OB 

DO  80  IV=IV1 ,IV2, IVT 

INITIALIZE  THE  NUfBER  OF  FALSE  ALARMS 

RMIS=0 

PERFORM  CLASSIFICATIONS  NEX  TIMES 
DO  888  IEX=1 , NEX 

IF (IDIS .  EQ.  1 ) PRINT  *, ' ************************************  1 

ADD  NOISE  TO  THE  ORIGINAL  DATA  FILES  TO  GENERATE 
TEST  TARGETS 

VM=10.**(FL0AT(IV)/10.) 

CALL  MN(VM, IS3 ,IS4,C, NC, NF, NS, IYA, IYW, IEX, SAP ) 

CALCULATE  NOISY  AMPLITUDE  AND  PHASE  RETURNS 

DO  921  1=1,  NS 
DO  921  J=1  ,NF 
AN  ( I ,  J  ]  =ABS  ( NC  ( I ,  J ) ) 

PN(I , J)=ATAN2(AIMAG(NC(I, J) ) , REAL(NC(I , J) ) ) 


CALCULATE  NOISY  W'S 


DO  922  1=1,  NS 
DO  922  J=1,NF-1 
PC=PN ( I , J ) 

DP=PN  ( I ,  J )  -PN ( I ,  J+1 ) 

I F { DP . GT . PI ) PC=PC-TWOPI 
I F ( DP . LT . -PI ) PC=PC+TWOPI 
WNd,  J)=WL(J)*PC-WL(J+1  )*PN(I,  J+1 ) 
IF(IRA.EQ.O)GO  TD  398 

CALCULATE  NOISY  RELATIVE  AMPLITUDES 

DO  349  IR1=1,NS 
DO  349  IR2=1,NF-1 

AN(IR1 ,IR2)=AN(IR1 ,IR2)/AN(IR1 ,IR2+1 ] 

COMPUTE  NEAREST-NEIGHBOR  DISTANCES 
BETWEEN  ANY  TWO  CLASSES 

CALL  DIST (IRA, IYW, IYA, A, W, AN, WN, NF, NS, RR,RK, D) 

PERFORM  NEAREST-NEIGHBOR  CLASSIFICATION 
AND  COMPUTE  THE  NUWER  OF  FALSE  ALARMS 

IF  ASPECT  ANGLES  ARE  ASSUMED  KNOWN  THEN 
CONSIDER  ONLY  THOSE  CLASSES  WHOSE  ASPECT  ANGLES 
ARE  THE  SAME  AS  THAT  OF  THE  TARGET,  OTHERWISE 
CONSIDER  ALL  CLASSES 

KASP= ( MAXASP-MI NASP ] /INCASP+1 
IFdAK.BQ.O)KASP=1 

FIND  THE  NEAREST-NEIGHBOR  OF  THE  TEST  TARTETS ' 

DO  7263  IKX=1  ,KASP 
DO  740  IT=IKX, NS, KASP 
DMIN=1 0000000000. 

DO  860  I=IKX, NS, KASP 
IF[D[I,IT).GE.DMIN)GO  TD  860 
DMIN=D(I,IT] 

IMIN=I 

CONTINUE 

IFdDIS.EQ.1  (PRINT  '  TARGET  :’,IT,’  MIN  =  ',1 
IF (IDIS.BQ.1  (PRINT  730,  (Ddl , IT )  ,11=1  ,NS) 
IFdDIS.  EQ.  2 ) WRITE (1 ,730)  CDCI1  , IT)  ,11=1 , NS) 
FORMAT (6X.15F5.0,/) 

IF  A  TARGET  IS  MISCLASSIFIED  INCREASE  THE 
NUMBER  OF  FALSE  ALARMS  BY  ONE 

I F (I MI N . NE . IT ) RMIS=RMIS  +1 

CONTINUE 

CONTINUE 

CALCULATE  THE  MAXIMUM  LIKELIHOOD  ESTIMATE 
OF  THE  ERROR  PROBABILITY 


RMIS=RMIS/ ( NEX*NS ) *1  DO 


PRINT  1113,1V, RMIS 

FORMAT ( 1  NOISE  POWER  =  ',13,';  X  ERRORS  =  \F5.1) 
CONTINUE 

TYPE  *,'00  YOU  WANT  TO  PROCEED  ?  Y=1  N=0' 

ACCEPT  *, IYN 
IFdYN.EQ.1 )G0  TO  555 
STOP 
ENO 

THIS  SUBROUTINE  CALCULATES  THE  NEAREST-NEIGHBOR 
DISTANCES  BETWEEN  ANY  TWO  CLASSES  IN  THE  DATA  BASE 

SUBROUTINE  DIST(IRA,IYW,IYA,A,W,AN,WN,NF,NS,RR,RK,D) 
REAL  D(1 8,1 8) ,  A(1 8,50)  ,W(  1  8,50] ,  AN(1  8, 50)  ,WN[1  8,50) 

DO  20  1=1  ,  NS 
DO  20  IT=1 ,NS 
D(I,IT)=0. 

DO  23  K=1 ,  NF-1 

IF(IYW.EQ.1)D(I,IT)=D(I,IT)+[ (W(I, K)-WN(IT,K) )*RK)**2 
IFdYA.EQ.I  )0(I , IT)=D(I , IT] +(A(I ,  K]-AN(IT, K]  )**2 
CONTINUE 

IF  ABSOLUTE  AMPLITUDES  ARE  USED  THEN  TAKE  INTO  ACCOUNT 
THE  LAST  FREQUENCY. 

IFdYA.EQ.I. AND.  IRA.  BQ.O)D(I,IT)=D(I,IT)  + 

1  (A[I,NF]-AN(IT,NF) )**2 
0(I,IT)=SQRT(0d,IT) ) 

CONTINUE 

RETURN 

END 

THIS  SUBROUTINE  ADDS  GAUSSIAN  NOISE  TO  THE  ORIGINAL 
DATA  FILES  TO  GENERATE  TEST  TARGETS 

SUBROUTINE  MN(Vfi,IS3,IS4,C,NC,NF,NS,IYA,IYW,IEX,SAP) 
COMPLEX  NC(1 8,50) ,C(1 8,50] 

REAL  SAP  (18) 

THE  GAUSSIAN  NOISE  IS  ADDITIVE  AND  HAS  VARIANCE  VR 
AND  MEAN  0. 

ASSUME  EQUAL  AMOUNTS  OF  NOISE  ARE  ADDED  TO 
THE  REAL  PARTS  AND  THE  IMAGINARY  PARTS 
OF  THE  DATA  POINTS 

VR1  =VR/2. 


DO  10  1=1 ,  NS 

RESET  THE  SEEDS  FOR  EACH  TEST  TARGET 


IS3=IS3-4 
IS4=IS4-6 
DO  10  J=1  ,NF 

GENERATE  TWO  GAUSSIAN  NUM3ERS 

CALL  GAUSS  ( IS3 , 0 . ,  VR1  ,  S3 ) 

CALL  GAUSS(IS4,0.,VR1,S4) 

GENERATE  A  TEST  TARGET 

NC(I,J]=CMPLX(REAL(C(I,J] ) +S3 , AIMAG(C(I , J) )+S4) 

IF  PHASE  RETURNS  ARE  UNKNOWN  THEN  SET  THEM  TO  ZERO 

IFdYA.ED.I  .AND.  IYW.  EQ.O)  NCd ,  J )  =CMPLX(ABS(C  (I ,  J )  )+S3,S4) 

RETURN 

END 

THIS  SUBROUTINE  READS  A  DATA  FILE  GENERATED  BY 
THE  B3.F0R  PROGRAM 

SUBROUTINE  DATA(MINASP,MAXASP, INCASP,FMIN,FMAX,FINC 
1,NF,NS,C,A,P,IRA] 

INTEGER *2  INFILE [ 5 ) ,FLNAM(5] .ASPI13) ,AIR[4) ,SHIP(5) 
COMPLEX  C(1 8,50] 

REAL  A(1  B,50]  ,P(1  8,50]  ,W(1  8,50] 

REALM  AM(  801 ) ,  PH  t  801  ] 

DATA  ASP/'00','15','30','45','6a'  ,'?5','90'  , 

1  ins'  <?n<  1  ^  s  '  i  cn  '  ICRI  i  on  i  / 

DATA  SHIP/' AD' , '  LR' , '  IJB  '  , '  KA'  ,  '7g’/,AIR/'B7  ' ,  'CO'  , '  C5 '  ,  't 
IBASE=0 

TYPE  AIRCRAFT  CLASSIFICATION  - >  TYPE  O' 

TYPE  *, '  SHIP  CLASSIFICATION  - >  TYPE  1 ' 

ACCEPT  *, IAS 

IF (IAS. EQ.O) PRINT  *, 'CLASSIFICATION  OF  AIRCRAFT' 

IFdAS.EQ.I )PRINT  *, 'CLASSIFICATION  OF  SHIPS’ 

IF(IAS.EQ.0)NST=4 

IFtIAS.EQ.1 )NST=5 

DO  772  1=1  ,NST 

IF (IAS. EQ. 01FLNAM1I )=AIR[I ) 

IF  (I  AS.  HI.  1 JFLNAM1I  ]=SHIP(I ) 

I F  ( I NCASP .  ED .  0 )  I  NCASP=1 
IB C= ( MAXASP-MI NASP ) /I NCASP 
DO  40  1=1, NST 

DO  111  IAP1=MINASP,MAXASP, INCASP 

IBASE=IBASE+1 

IASP=IAP1/15+1 

INFILEd  )  =  'R: ' 

INFILE(2)=FLNAM(I ) 

INFILE(3)='FV' 

INFILE(4]=ASP(IASP) 

INFILE(5)=0 

CALL  TR (INFILE, AM, PH, TFMIN.TFINC.SF) 


C  NOTE  THE  DATA  FILES 

C  ARE  ALREADY  SCALED  BY  THE  SCALE  FACTOR 

C 

I MAX= { FMAX-TFMI NJ/TFINC+1 
IMIN= ( FMIN-TFMI N ) /TFINC+1 
IC=FINC/TFINC 
NF1=0 

DO  2111  J=IMIN,IMAX,IC 
C 

C  CONVERT  AMPLITUDE  SCALE  FORM  CM  TO  M 

C  BY  DIVISION  OF  100 

C 

NF1  =NF1  +1 

A  ( IB  ASE ,  NF1 )  =1 0 .  **  ( AM  l  J )  /20 .1/100. 

2111  PlIBASE,  NF1  )=PH  (J } 

111  CONTINUE 

40  CONTINUE 

NS=IBASE 

6664  DO  1436  1=1  ,  NS 
DO  1436  J=1,NF 

1436  Cd,J)=CMPLX(Ad,J)*COS(Pd,J) ),A[I,J)*SIN[P(I,J))] 
IFINF.LE.NF1  )G0  TO  1435 

DO  1437  1=1, NS 
Pd,NFJ=0 
All , NF)=0 

1437  C(I,NFJ=0 

1435  RETURN 

END 

C 

C  THIS  SUBROUTINE  CALCULATES  THE  VARIANCE  OF 

C  A  SET  OF  NUIGERS 

C 

SUBROUTINE  VARN(A,NF,NS,VA,INC) 

REAL  All  8,50) 

RMEAN=0 
DO  10  1=1, NS 
DO  10  J=1  ,  NF.INC 
10  RMEAN=A ( I , J ) +RMEAN 

RMEAN=RMEAN/ ( NS* ( NF/I NC ) ) 

VA=0 

DO  20  1=1,  NS 
DO  20  J=1, NF.INC 
20  VA=VA+fAd,  J)-RMEAN)**2 

VA=VA/(NS*(NF/INC) ) 

RETURN 

END 

C 

C  THIS  SUBROUTINE  GENERATES  A  GAUSSIAN  RANDOM  NUMBER 

C 

SUBROUTINE  GAUSS  dS,  AM,  V,SJ 
A=0.0 

DO  50  1=1,12 
A=A+RAN(IS) 

50  CONTINUE 


c 

C  THIS  SUBROUTINE  READS  A  FILE 

C 

SUBROUTINE  TR (INFILE, AM, PH , FMIN, FINC, SF) 

INTEGER*2  INFILE15) ,LINE2(24) 

INTEGER*4  ID(12) 

REAL*4  AM  ( 801 ) ,  PH  ( 801 ) 

INCLUDE  1 SYSSLIBRARY : FORIOSDEF ' 

8106  OPEN (UNIT=B,NAME=INFILE,TYPE=  'OLD' , IOSTAT=IERR, ERR=8100) 

READ ( 8,30) (ID(I), 1=1,1 2), LINE2.NP, FMIN, FINC, SF 
30  FORMAT ( 1 2A4 , 24A2 , 4A4) 

READ  [  8, 40 )( AM  (I) ,  1=1 ,  NP ),( PH  (I) ,  1=1 ,  NP) 

40  F0RMAT(128A4] 

GO  TO  331 

8100  IF(IERR.BQ.  FOR*IOS_FILNOTFOU)THEN 

TYPE  1112, INFILE 

1112  FORMAT ( 1  FILE  :  \15A2,'  WAS  NOT  FOUND',//, 

1'  ENTER  FILENAME  AGAIN') 

ELSE  IF  (IERR.EQ. FOR$IOS_FILNAMSPE)THEN 
TYPE  1113, INFILE 

1113  FORMAT ('  FILE  :  ',15A2,'  WAS  BAD,  ENTER  NEW  FILENAME') 

ELSE 

TYPE  *, 'UNRECOVERABLE  ERROR,  CODE=',IERR 

STOP 

ENDIF 

ACCEPT  1114, (INFILE (I IT) ,IIT=1 ,14) 

1114  FORMAT  (14A2) 

GO  TO  8106 

331  CLOSE (UNIT=8,DISP=' SAVE ' ) 

RETURN 

END 

C 

C  THIS  SUBROUTINE  DECODES  PARAMETERS  READ 

C  FROM  THE  DATA  FILES,  WHICH  ARE  CODED  IN  ASCII  CODES 

C 

SUBROUTINE  DCDEtNP, FMIN, FINC) 

COMMON  BUFF 
BYTE  BUFF(6656) 

INTEGERS  I  MIN,  IINC,  NP 
C 

C  NO  OF  DATA  POINTS  IS  STORED  IN  THREE  CHARACTERS,  AND 

C  STARTING  FREQ  AND  FREQ  INC.  IN  5  CHARACTERS 

C 

CHARACTER  *3  CNL 
CHARACTER *5  CFF.CINC 

EQUIVALENCE  (BUFF(124) ,CNLJ , (BUFFC131 ) ,CFF) , (BUFF1140) ,CINC) 
C 
C 

C  CONVERT  CHARACTERS  INTO  THEIR  NUMERICAL  EQUIVALENTS 

C 

DECODE ( 3 , 1 00 , CNL ) NP 
100  FORMAT (15) 

DEC00E(5,100,CFF)IMIN 
DECODE (5, 100, Cl NC) I INC 
FMIN=FL0AT(IMIN)/1000. 
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c 

C  PROGRAM  NAME  :  TIME. FOR 

C 

C 

C  THIS  PROGRAM  CLASSIFIES  TARGETS  IN  THE  TIME  DOMAIN 

C 

C 

COMPLEX  C  ( 1  8, 50 ) ,  NC  ( 1  8, 50 ) 

COMPLEX  FFFT { 32 ] 

REAL  D(1 8,1 8)  ,SFFT ( 8) 

REAL  SP CIS)  ,SNP(1 8} 

COMMON  /B  LK1 /SFFT , FFFT 

TYPE  *,'  ••***  YOU  MUST  ASSIGN  DATA  DIRECTORY  TO  R  **•**’ 

TYPE  WANT  TO  PRINT  DISTANCES  IN  CLASSIFICATION  ?  Y=1  N=0' 
555  ACCEPT  *,IDIS 

C 

C  INPUT  TWO  ARBITRARY  LARGE  ODD  INTEGERS  WHICH 

C  WILL  SERVE  AS  SEEOS  FOR  GENERATING  RANDOM  NU  USERS 

C 

TYPE  ‘.'INPUT  2  ARBITRARY  LARGE  ODD  INTEGERS' 

ACCEPT  *, IS3.IS4 
C 

C  SPECIFY  THE  MINIMUM,  MAXIHJM,  AND  INCREMENT  OF  THE 

C  ADDITIVE  GAUSSIAN  NOISE  POWERS 

C 

TYPE  *, 'INPUT  MIN, MAX, INC  OF  THE  NOISE  POWER  [DB)> 

ACCEPT  *,IV1 .IV2.IVT 
C 

C  SPECIFY  TARGET  ASPECT  ANGLES 

C 

TYPE  *, 'INPUT  MIN, MAX, INC  OF  ASP' 

ACCEPT  *,MINASP,MAXASP,INCASP 
C 

C  SPECIFY  FREQUENCIES 

C 

TYPE  *,' INPUT  FMIN.FINC.NO  OF  FREQ' 

ACCEPT  *, FMIN, FINC, NF 
C 

C  SET  THE  NUffiER  OF  CLASSIFICATION  EXPERIMENTS 

C 

TYPE  *, 'INPUT  NO  OF  EXPERIMENTS' 

ACCEPT  *, NEX 
C 

TYPE  *,  'ASP  KNOWN  — >  TYPE  1' 

TYPE  *, 'ASP  UNKNOWN  — >  TYPE  O' 

ACCEPT  *,  IAK 

FMAX=FMIN+(NF-1 )*FINC 

PRINT  V  LARGE  INTEGERS  :  ',IS3,IS4 

IF ( IAK. BQ.1) PRINT  *, '  ASP  KNOWN' 

IF ( IAK. EQ.O) PRINT  ♦,'  ASP  UNKNOWN' 

PRINT  V  ADDITIVE  NOISE  ADDED1 

PRINT  V  CLASSIFICATION  IN  TIME  DOMAIN' 

PRINT  *,'  NO  OF  EXPERIMENTS  ',NEX 
PRINT  2B39,FMIN, FMAX, FINC, NF 

2839  FORMAT l  ’  FMIN  =  ',FB.3,'  FMAX  =  \FB.3,'  FINC  =  ', 
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V*  '  j 


sVvAV'*  »■ 


.  W  -  »  *•  •  <1 


1F8.3,/,’  NO  OF  FREQ  =  ',13] 

PRINT  2849,MINASP,MAXASP,INCASP 
2849  FORMATC  MIN, MAX, AND  INC  OF  THE  ASPECTS  :  ’,315) 

C 

C  READ  DATA  FILES 

C 

CALL  QATA{MINASP,MAXASP,INCASP,FMIN,FMAX,FINC, 

INF, NS, C, IRA) 

C 

C  CALCULATE  A  COSINE  TABLE  FOR  FAST  FOURIER  TRANSFORMS 

C  BY  CALLING  A  SUBROUTINE  IN  THE  SCIENTIFIC  LIBRARY 

C 

CALL  F0RT(FFFT,5,SFFT,0,IFERR) 

C 

C  CALCULATE  SQUARE-ROOTS  OF  TOTAL  SIGNAL  POWERS 

C 

DO  9B274  1=1 ,  NS 
SP(I)=0 

DO  9824  J=1 ,  NF 

9B24  SP(I  )=SP(I )  +ABS1C  (I ,  J )  ]**2 
9827  4  SP(I)=SQRT(SPtin 

C 

PRINT  *, 1  NO  OF  TARGETS  ',NS 
PRINT  18220 

18220  FORMAT (IX  I*************************************** 1 ) 

C 

C  VARY  THE  NOISE  POWERS  FROM  IV1  TO  IV2  DBS 

C  AT  INCREMENT  OF  IVT  DBS 

C 

DO  80  IV=IV1 ,IV2,IVT 
C 

C  INITIALIZE  THE  TOTAL  MJ43ER  OF  FALSE  ALARMS 

C 

RMIS=0 

C 

C 

C  PERFORM  CLASSIFICATION  EXPERIMENTS  NEX  TIMES 

C 

DO  888  IEX=1 ,NEX 

IFfIDIS.EQ.1 ) PRINT  •,•••••***************•********** 

C 

C  ADO  NOISES  TO  ORIGINAL  DATA  FILES  TO  GENERATE 

C  TEST  TARGETS 

C 

VM=10  .•*(FL0AT(IV)/10. ) 

CALL  MN(VM, IEX, IS3  TS4,C, NC, NF, NS, SP) 

C 

C 

C  CALCULATE  SQUARE-ROOTS  OF  TOTAL  POWERS  OF 

C  THE  NOISY  TARGETS 

C 

00  9221  1=1,  NS 
SNPII )=0 
DO  921  J=1  ,NF 

SNPd  )=SNP(I )  +ABS(  NC  [I ,  J )  )**2 


921  CONTINUE 

9221  SNP(I)=SOFrr[SNP(I)) 

C 

C  COMPUTE  CROSS-COEFFICIENTS  BETWEEN  ANY  TWO  CLASSES 

C 

fB=INT(FMIN/FINC) 

479  CALL  DIST (C, NC,SP,SNP, NF, NS, D, (B ] 

C 

C  PERFORM  TlfC-DOMAIN  CLASSIFICATIONS, 

C  AND  COMPUTE  NO  OF  FALSE  ALARMS 

C 

C  IF  ASPECT  ANGLES  ARE  ASSUMED  KNOWN  THEN 

C  CONSIDER  ONLY  THOSE  CLASSES  WHOSE  ASPECT  ANGLES 

C  ARE  THE  SAME  AS  THAT  OF  A  TARGET,  OHTERWISE 

C  CONSIDER  ALL  CLASSES 

C 

KASP=  ( MAXASP-MI NASP1  /INCASP+1 
IF (IAK. EQ. 0) KASP=1 
C 

00  7  263  IKX=1  ,KASP 
DO  740  IT=IKX,NS,KASP 
DMAX=-1 0000000000. 

C 

C  FIND  A  CLASS  WICH  HAS  THE  MAXIHJM  CORRELATION 

C  WITH  A  TARGET 

C 

DO  060  I=IKX,NS,KASP 

IF (D (I, IT) . LE.DMAXJGO  TO  860 

DMAX=0 ( I , IT ) 

IMAX=I 

860  CONTINUE 

C 

IFdDIS.BQ.DPRINT  *, '  TARGET  :',IT,'  MAX  =  \IMAX 
IFIIDIS.BQ.I ] PRINT  730,(0(11 , IT) ,11=1 , NS) 
IF(IDIS.EQ.2)WRITE(1 ,730) (0(12, IT) ,12=1 , NS) 

730  FORMAT (6X.15F5.2,/) 

C 

C  IF  A  TARGET  IS  MISCLASSIFIED  THEN  INCREASE  THE 

C  WJfCER  OF  FALSE  ALARMS  8Y  1 

C 

740  IF(IMAX.NE.IT)RMIS=RMIS+1 

C 

7263  CONTINUE 

880  CONTINUE 

C 

C  CALCULATE  THE  MAXI  HIM  LIKELIHOOO  ESTIMATE  OF 

C  THE  ERROR  PROBABILITY 

C 

RMIS=RMI S/FLOAT ( NEX*NS ) *1 00 . 

C 

C 

PRINT  1113,1V, RMIS 

1113  FORMAT  ( '  NOISE  POWER  =  ',15,';  *  ERRORS  =  \F5.1) 

80  CONTINUE 

TYPE  *, 'DO  YOU  WANT  TO  PROCEED  ?  Y=1  N=0’ 
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ACCEPT  *, IYN 
IFUYN.EQ.  1  )G0  TO  555 
STOP 
END 

THIS  SUBROUTINE  CALCULATES  THE  CROSS-CORRELATION 
COEFFICIENTS  BETWEEN  TWO  CLASSES 

SUBROUTINE  DIST (C, NC, SP, SNP, NF,  NS, D,  ] 

REAL  D(1 8,1  8]  ,SP(1 8) , SNP  11 8)  ,S( 8) 

COMPLEX  C(1 8,50) ,NC(1 8,50) ,RF(50) ,F(32) 

DO  20  1=1,  NS 
DO  20  IT=1 ,NS 
DO  23  K=1  ,NF 

CALCULATE  PRODUCTS  OF  RADAR  RETURNS  OF  A  CLASS 
AND  COMPLEX  CONJUGATES  OF  A  TARGET'S  RETURNS 

RF[K)=C(I,K)*CONJG(NC(IT,K) 1 

CONTINUE 

CALCULATE  THE  INVERSE  FAST-FOURIER 
TRANSFORM  (  IFFT  )  OF  THE  PRODUCTS 

CALL  IDIFTIRF, NF, fe , STM) 

NORMALIZE  THE  OUTPUTS  OF  THE  IFFT 

0(I,IT)=STM/(2.*SP(I)*SNP(ITJ) 

CONTINUE 

RETURN 

END 

THIS  SUBROUTINE  ADOS  GAUSSIAN  NOISE  TO  THE  ORIGINAL 
DATA  FILES  TO  GENERATE  TEST  TARGETS 

SUBROUTINE  MN(VR,IEX, IS3 , IS4,C, NC, NF, NS, SP) 

COMPLEX  NCC1 8,50) ,C(1 8,50) 

REAL  SP(1 8) 

THE  GAUSSIAN  NOISE  ARE  ADDITIVE  AND  HAS  VARIANCE  VR 
AND  MEAN  0. 

VR1=VR/2 
00  10  1=1, NS 

RESET  SEEDS  FOR  EACH  CLASS 

IS3=IS3-4 
IS4=IS4-2 
DO  10  J=1  ,  NF 

GENERATE  TWO  GAUSSIAN  NU)CERS 
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CALL  GAUSS ( IS3 , 0 . , VR1 , S3 ) 

CALL  GAUSS{IS4,0. ,VR1 ,S4) 

C 

C  GENERATE  A  TEST  TARGET'S  RETURN 

C 

10  NC(I,J)=CMPLX(REAL(C(I,J) ] +S3 , AIMAGfCtI , J)) +64) 

C 

RETURN 

END 

C 

C  THIS  SUBROUTINE  READS  DATA  FILES  GENERATED  BY  PROGRAM 

C  83. FOR 

C 

SUBROUTINE  DATA(MINASP,MAXASP, INCASP, FMIN, FMAX.FINC 
1 ,NF,NS,C,IRA) 

INTEGER*2  INFILE (5) ,FLNAM(5) ,ASP(13) ,AIR[4) .SHIP [5) 

COMPLEX  C(1 8,50) 

REAL  A(18,50)  ,P(18,50)  ,W(1B,50) 

REALM  AMI  801 )  , PH  (801  ) 

DATA  ASP/' 00' ,  '15' ,  '30' ,  *45' ,  '60' ,  *75' ,  ’  90' , 

1  '05'  ,  '20'  ,  '35' ,  '50'  ,  '65' ,  '  BO'/ 

DATA  SHIP/'AD'  ,  'LR'  ,  'LB'  ,  'KA'  ,  '  79'/,  AIR/ 'B7  ' , '  CO '  ,  '  C5  '  ,  'DC'/ 
IBASE=0 

TYPE  *,'  AIRCRAFT  CLASSIFICATION  - >  TYPE  O' 

TYPE  SHIP  CLASSIFICATION  - >  TYPE  1' 

ACCEPT  *, IAS 
C 

C  NST  =  NO  OF  TARGETS 

C 

IF (IAS . EQ. 0) NST=4 
IF ( IAS . ED. 1 ) NST=5 
C 

C  GENERATE  FILE  NAMES 

C 

DO  772  1=1,  NST 

IF (IAS. EQ. 0 ) FLNAMt I ) =AIR ( I ) 

772  IF11AS.EQ.1)FLNAM(I)=SHIPII) 

IF  (INCASP.  HI.  0)  INCAS P=1 
DO  40  1=1  ,  NST 

00  111  IAP1=MINASP,MAXASP, INCASP 

IBASE=IBASE+1 

IASP=IAP1/15+1 

INFILEtl )  =  'R:  ' 

INFILE(2!=FLNAM(I) 

INFILE(3)='FV’ 

INFILE(4)=ASP(IASP) 

INFILE(5)=0 

C 

C  READ  A  FILE 

C 

CALL  TR (INFILE, AM, PH, TFMIN,TFINC,SF) 

C 

C  NOTE  THE  INPUT  FILES 

C  HAVE  ALREADY  BEEN  SCALED  BY  THE  SCALE  FACTOR 
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I MAX= ( FMAX-TFMI N ) /TFINC+1 
I MI N= ( FMI N-TFMI NJ/TFINC+1 
IC=FINCAFINC 
NF1=0 

00  2111  J=IMIN,  IMAX,  IC 
C 

C  CONVERT  AMPLITUDE  SCALE  FORM  CM  TO  M 

C  BY  DIVISION  OF  100 

C 

NF1=NF1  +1 

AtIBASE,  NF1 )  =10.** (AM (J 1/20.  ]/100. 

2111  PIIBASE,NF1)=PH(J) 

111  CONTINUE 

40  CONTINUE 

NS=IBASE 

6664  DO  1436  1=1 , NS 
DO  1436  J=1  ,NF 

1436  C ( I , J)=CMPLX(A(I , J) *C0StP(I,  J) ) , A(I,J) *SIN{P(I , J) ) } 
IFtNF.LE.NF1  )G0  TO  1435 

DO  1437  1=1, NS 
PlI , NF)=0 
A(I,NF)=0 

1437  C1I ,  NF)=0 

1435  RETURN 

END 

C 

C  THIS  SUBROUTINE  GENERATES  GAUSSIAN  RANDOM  NUMBERS 

C  WHICH  HAVE  A  VARIANVE  V  AND  MEAN  AM 

C 

SUBROUTINE  GAUSSIIS, AM, V,S) 

A=0.0 

DO  50  1=1  ,12 
A=A+RAN(IS) 

50  CONTINUE 

S=(A-6 . )*SQRT(V) +AM 
RETURN 
END 
C 

C  THIS  SUBROUTINE  CALCULATES  THE  INVERSE 

C  FOURIER  TRANSFORM  OF  A  SET  OF  COMPLEX  DATA  POINTS 

C  ANO  FIND  THE  MAXIHJM  OF  THE  RESULTS 

C  IT  ASSUMES  THAT  THE  DATA  POINTS  ARE  SYMMETRIC 

C  AROUND  ZERO  FREQUENCY 

C 

SUBROUTINE  IOIFTtRF, NF,  fB , STM ) 

COMPLEX  F(32) ,RF(50) 

REAL  St 6) 

COMMON  /BLK1/S.F 
NP=32 

DO  40  1=1,  NP 
40  FlI 1=10- ,0.1 

K=0 

DO  50  I=fB+1,fB+NF 
K=K+1 
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